Theoretical Model and Numerical Simulation of Adsorption and

Article Cite This: J. Phys. Chem. C 2018, 122, 9465−9473

pubs.acs.org/JPCC

Theoretical Model and Numerical Simulation of Adsorption and Deformation in Flexible Metal−Organic Frameworks Sahar Bakhshian and Muhammad Sahimi*

Downloaded via UNIV OF LOUISIANA AT LAFAYETTE on September 28, 2018 at 14:51:22 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

Mork Family Department of Chemical Engineering and Materials Science, University of Southern California, Los Angeles, California 90089-1211, United States ABSTRACT: Metal−organic frameworks (MOFs) have recently attracted considerable attention as new nanoporous materials with applications to separation of fluid mixtures, catalysis, gas capture and storage, and drug delivery. A subclass of MOFs, the MIL-53 (with Al or Cr) family, exhibits a complex structural phase transition, often referred to as the breathing transition, which occurs when gas molecules are adsorbed in its pore space. In this phenomenon, the material morphology oscillates between two distinct phases, usually referred to as the narrow-pore (np) and large-pore (lp) phases. We describe a statistical mechanical model based on the energetics of the system that couples the adsorbates with the host and the ensuing structural deformation of the materials. The Hamiltonian of the system consists of the elastic energy of the MOF, the fluid phase energy, and the interaction energy between the gas and MOFs. Minimizing the total energy with respect to both the gas density and the displacement field in the material yields the governing equations for both. The model is used to study the adsorption of CO2 and CH4 in MIL-53(Al). The results demonstrate that all reported experimental features of the phenomenon are produced by the model. In particular, consistent with experimental data, the model predicts that the material experiences contraction in the np phase, whereas it undergoes swelling in the lp phase, and that the two phases coexist in the transition zone. The model is, however, completely general and is applicable to any type of guest-responsive material that undergoes deformation as a result of its interaction with the guest molecules.

1. INTRODUCTION Metal−organic frameworks (MOFs), also known as porous coordination networks, constitute a class of porous crystalline materials that are formed by metallic species covalently linked together through organic linkers.1−3 So far, a large range of MOF materials have been synthesized through various combinations of metal ions and organic elements.4−7 However, because of the combinatorial possibilities, there are millions of potential MOF materials. Thus, large-scale molecular simulations have been used8−11 in order to evaluate their various properties, such as porosity and pore surface areas, and assess their potential applications to such important problems as hydrogen storage and CO2 adsorption. The type of metal and organic subunits greatly influences the properties of MOFs. Such materials have several important characteristics, including large surface area, tunable porosity, and the ability to undergo reversible deformation, which explains why they are finding a remarkable number of applications in gas storage, separation of fluid mixtures, and adsorption technologies.12−19 In particular, some MOFs exhibit structural transition and deformation with large amplitudes when they are exposed to such external stimuli as light, temperature, electric field, and pressure or take part in gas adsorption or liquid absorption.4,20,21 Among them is a subset of MOFs called MIL-53 (where MIL stands for Materials of Institüt Lavoisier). They consist of corner metals, such as Al, © 2018 American Chemical Society

Cr, Ga, and Fe, which are connected with benzenedicarboxylate ligands, giving rise to essentially one-dimensional diamondshaped pore channels.2,22,23 They represent flexible MOFs, and when exposed to gas adsorbates, they undergo what is usually referred to as breathing phenomenon,4 which is a reversible flexing or deformation of the framework as a function of the amount and type of the adsorbates. In particular, the MIL-53 structure undergoes a transformation from a large pore (lp) to a narrow pore (np) phase and switches back to the lp state during the breathing phenomenon when exposed to such adsorbate molecules as CO2 and CH4. This implies that multiple stable states can exist in a single framework. Such unusual and useful characteristics of this class of MOFs have provided the impetus for researchers to consider innovative applications for flexible MOFs. Several studies have demonstrated step-wise adsorption isotherms in some MOFs, including MIL-53 and MIL-88. Such a complex behavior is attributed to the occurrence of a sort of a first-order phase transition (see below), including breathing, or crystal-to-crystal transformation in the structure of MOFs. In particular, it has been shown that MIL-88 may swell up to 200% in cell volume during fluid adsorption, whereas the Received: January 27, 2018 Revised: March 14, 2018 Published: April 9, 2018 9465

DOI: 10.1021/acs.jpcc.8b00924 J. Phys. Chem. C 2018, 122, 9465−9473

Article

The Journal of Physical Chemistry C amplitude of breathing may reach about 38% in MIL-53, when it adsorbs CO2.24−26 It is clear that such phenomena are due to the interaction of the guest molecules with the structure of MOFs. The aforementioned deformation and structural transition in flexible MOFs that are induced by adsorption have been studied both experimentally and theoretically.25−35 Coudert et al.36 studied gas adsorption and the structural transition that it induces in such materials by formulating the problem based on the thermodynamics of the system and calculated the variation of the free energy of the system during the phase transition in the structure of the material. Neimark et al.37 proposed a model of the breathing transition in MIL-53 based on the variation of the stress during adsorption of a gas. Their model predicted coexistence of the aforementioned np and lp phases at pressures close to the one at which the breathing transition occurs. Neimark et al. also studied22 adsorption-induced structural transition in flexible MIL-53(Cr). They found that adsorption of guest molecules in such materials induces small reversible elastic deformation of around 2 and 4% in, respectively, the lp and np phases. In addition, their model predicted an irreversible plastic deformation of about 40% in the structure during its transition from the np to lp phases. Another model for adsorption-induced breathing transition in MIL-53(Cr) was proposed by Ghysels et al.4 Using a parameterized energy model, they derived analytical expressions for the free energy of the guest molecules, the host material, and the host−guest interaction as a function of cell parameters. In this paper, we describe a statistical mechanical model of adsorption-induced deformation in flexible MOFs. The model is based on a coupling between an elastic porous solid and its interaction with a gas that resides and is adsorbed in the material pores. A general expression for the Hamiltonian of the system is used that takes into account not only the interaction of the gas molecules with the material but also the interaction between the guest molecules themselves. The model is then used to study adsorption of CO2 and CH4 in MIL-53(Al) at various temperatures and the structural transitions that occur in it. The rest of this paper is organized as follows. The theoretical model of adsorption and deformation is described in the next section. Next, the procedure for the numerical simulation of the theoretical model is described, after which the results are presented and discussed. The paper is summarized in the last section.

Hf = kBT ∑ [ρi ln ρi + (1 − ρi )ln(1 − ρi )] i

1 − εff ∑ ∑ ρi ρi + Z − μ b ∑ ρi 2 i Z i

where ρi is the average fluid occupancy (or average density) at lattice site i, T is the temperature, kB is the Boltzmann constant, Z is the lattice coordination number, μb is the bulk fluid chemical potential, and ϵff is a parameter that describes the interaction between the fluid molecules. We make the typical assumption that the interactions are only between the nearest neighbor molecules in the lattice. To compute the fluid occupancy and, hence, its density, we define for each lattice site (m) i two indicator functions such that I(f) =1 i = 0 and 1 and Ii and 0 represent, respectively, the fluid and the material. Then, (m) the fluid density at each site is given by ρi = ⟨I(f) i [1 − Ii ]⟩. The expression for the energy Hm depends, of course, on the material and its constitutive equation. We assume that the MOF is linearly elastic. The assumption can, of course, be relaxed if a more precise constitutive equation for MOFs becomes available. For a linearly elastic material, the energy Hm is given by the standard expression Hm =

1 ∑ (K i∇·ui 2 + 2σi∇·ui) 2 i

(3)

where Ki is the elastic stiffness. Here, ∇·ui is the strain in the element i with ui being the displacement at i and σi is the corresponding stress. It remains to specify the energy Hfm that represents the strength of the interaction between the fluid and the material and couples the two. Two factors contribute to Hfm. One is the wetting energy Hw at the fluid−solid interface, which is clearly proportional to the total fluid density. Thus, Hw = − 1 ϵfm∑iρi, 2 where ϵfm represents the strength of the interaction between the fluid and the material. The second contributing factor is Hd, the change in the energy of the material as a result of deformation induced by adsorption of the gases. Hd is given by41,42 Hd =

1 λ∑ ρ ∇·u i 2 i i

(4)

with λ being the strain energy parameter that couples the change in the strain energy of the material as a result of its deformation, induced by adsorption of the fluid. Note that the sum runs over all sites of the computational grid. Equation 4 completes the formulation of the problem. Let us point out that formulation of the problem of gas sorption in porous materials by energy minimization has been done in the past. What sets apart the model is its generality and flexibility. For example, one does not have to assume that the material is linearly elastic. Thus, if, for example, we wish to apply the model to liquid absorption by polymers and the resulting deformation, we can use a suitable expression for the total energy of the polymer, which does not have to be linearly elastic and represents better its deformation. In addition, for absorption of liquids by polymers and similar materials, the lattice gas model would not be suitable, and one must use another equation of state. The model makes no assumption regarding either aspect, and thus it is quite general. Other existing models make several assumptions in order to arrive at a practicable model. In addition, the model is applicable to any

2. MODEL DESCRIPTION We assume that the deformation of MOFs that is induced by adsorption of CO2 or CH4 is a slow, quasi-static phenomenon. Then, the overall Hamiltonian of the systemthe total energy H of the porous material and its fluid contentis given by H = Hf + Hm + Hfm

(2)

(1)

Here, Hf, Hm, and Hfm represent, respectively, the energies that are attributed to the fluid, the solid material, and the interactions between the two. Because the fluid that we consider in this paper is either CO2 or CH4, which are relatively simple molecules, the energy associated with them is expressed by a mean-field density functional theory (DFT):38−40 9466


Article

The Journal of Physical Chemistry C type of microstructure provided that its morphology is wellcharacterized and known, or its image is available, so that the computations can be carried out directly with the image.

1 2

∑ K nl ul = ∑ eniλ ρ l

i

l

(9)

which is the analogue of eq 6 for the entire computational grid. Equations 5 and 9 are then solved numerically in order to determine the displacement and fluid density field. Setting up the global stiffness matrix K is a standard part of any FE computation. One needs the connectivity of the grid elements, the number of degree of freedoms per node and nodes per element, the number of elements, the spatial distribution of the elements, and the element stiffness matrices. The connectivity is a set of data associated with an element that determines which nodes (or vertices) form the element. This is typically formed by moving about the element counterclockwise and registering the node numbers as they are met. The crystalline structures of the np and lp phases with the cell parameters a, b, and c are shown in Figure 1. The view is

3. COMPUTER SIMULATION The gaseous adsorbate is assumed to be in contact with a bulk reservoir at chemical potential μb and temperature T. When the fluid inside the pores reaches equilibrium with MOF, the Hamiltonian is minimum. Thus, to obtain the governing equations for the fluid density and strain ui throughout the system, we minimize H, that is, we set, ∂H/∂ρi = ∂H/∂ui = 0. The two sets of coupled equations are solved numerically, for which we use the finite-element (FE) method (see below). Minimizing the Hamiltonian H with respect to ρi and ui yields the following equation:41,42 −1 ⎧ ⎡ ⎛ ⎞⎤⎫ ⎪ ⎪ ⎥ ⎢ ρi = ⎨1 + exp −β ⎜⎜μ b + εfm + εff ∑ ρi + Z − λ∇·u i⎟⎟ ⎬ ⎪ ⎢⎣ ⎝ ⎠⎥⎦⎪ ⎭ ⎩ Z

(5) −1

for the fluid density ρi, where β = (kBT) , and

k eul =

1 eniλρi 2

(6)

for the displacement ul of vertex l in element i of the computational grid. The unit vectors eni point from the center of element i to its vertices n. ke, the elemental stiffness matrix, is given by ke =

∫V UTCU dV

Figure 1. Schematic of the MIL-53 framework with the cell parameters a, b, and c. The view is along the c axis. (a) lp and (b) np structures.

(7)

along the c axis [perpendicular to the (a,b) plane]. For the lp and np shapes, the parameters used were,4 respectively, a/b = 1.3 and α = 14° and a/b = 2.7 and α = 46°. The shear modulus G and the Poisson’s ratio ν were taken to be, respectively, 7.5 GPa and 0.2 for the np phase and 1.5 GPa and 0.2 for the lp phase, which are the typical values reported22 for the MIL-53 MOF. Figure 2 presents the 3D images of the MOF structure that were discretized with adaptive tetrahedral elements. The

Here, C and U denote, respectively, the elastic coefficient and the strain-displacement matrices and superscript T indicates the transpose operation. The integral in eq 7 is over the volume V of the elements of the computational grid. C is given by ⎛ χ + 2μ ⎜ ⎜χ ⎜χ C=⎜ ⎜0 ⎜ ⎜0 ⎜ ⎝0

0 0 0⎞ ⎟ χ + 2μ χ 0 0 0⎟ χ χ + 2μ 0 0 0 ⎟ ⎟ μ 0 0⎟ 0 0 ⎟ 0 0 0 μ 0⎟ ⎟ 0 0 0 0 μ⎠ χ

χ

(8)

The two Lamé constants χ and μ are related to the shear modulus G and the Poisson’s ratio ν of the material by χ = 2νG/(1 − 2ν) and μ = G. The MOF structure is discretized using tetrahedron grid blocks. The normalized coordinates, (ξ,η,ζ) ≡ (f1,f 2,f 3), which are related to the Cartesian coordinates (x, y, z) by x = x1 + ∑4i=2(xi − x1)f i−1, y = y1 + ∑4i=2(yi − y1)f i−1, and z = z1 + ∑4i=2(zi − z1)f i−1, are introduced with (xi, yi, zi) (i = 1−4) representing the coordinates of the tetrahedron vertices. Using a shape vector N = (N1,N2,N3,N4)T with its components Ni(ξ,η,ζ) given by N1 = 1 − Σj = 24f j−1 and Nj = f j−1 with j = 2−4, the 3 × 4 matrix U in eq 7 is given by U = DN, with D being the differentiation operator matrix, so that the entries Uij are given by Uij = ∂Ni/∂f j. We then compute the global stiffness matrix K, given the stiffness matrices of the tetrahedral elements ke. The governing equations for the displacements of the elemental vertices are the following system of linear equations:

Figure 2. Discretized structures of the (a) lp and (b) np frameworks using tetrahedral cells.

computational grid must be resolved enough to represent accurately the morphology of the material and generate, with reasonable computational time, accurate adsorption isotherms and the resulting deformation. A computational domain with 478, 710 tetrahedra was used for simulating sorption and deformation for the np phase, whereas the corresponding grid 9467


Article

The Journal of Physical Chemistry C

Table 1. Dependence of the Parameters λ* and ϵ* on Temperature for CO2 in the np and lp Phases

for the lp phase consisted of 438 940 tetrahedral elements. Both were selected based on preliminary simulations, which indicated that the adsorption isotherms would not change with higher resolution grids with a larger number of smaller elements.

4. RESULTS AND DISCUSSION We have applied the model to sorption of CO2 and CH4 in flexible MOF-53(Al). In order to carry out the computations, however, we must first supply the parameters of the model. Thus, we first describe their estimation. 4.1. Parameter Estimation. The three parameters of the model, namely, the fluid−fluid interaction parameter ϵff, the fluid−material interaction parameter ϵfm, and the coupling strain energy parameter λ, may be estimated by various methods. For example, all three may be estimated by direct molecular dynamics (MD) simulations. Alternatively and more simply, ϵff and ϵfm may be taken to be the Lennard-Jones energy parameters for the interactions between, respectively, the fluid molecules themselves and between the fluid and the surface of the material. In that case, the parameter λ may also be estimated by MD simulation using the following method.42 Suppose that L is the length of a pore of the material. Knowing the chemical composition of the material, we first generate an atomistic model of the pore by the standard MD simulation, that is, by starting from an initial model of the pore and the chemical composition of its surface, minimizing the energy and carrying out MD simulation to obtain an equilibrium model. The gaseous adsorbate or the absorbing liquid is then introduced into the pore. MD simulations are then carried out in order for the pore and its content to reach equilibrium with the bulk state at a given chemical potential or, equivalently, pressure of the gas or liquid reservoir. The fluid molecules interact with the pore walls. Similar to an actual experiment, the chemical potential μ (the external pressure) of the bulk state is increased, which changes the distribution of the fluid in the pore as it tries to deform the pore length L. The change in L as a function of μ is computed by the MD simulation. According to eq 6, λ = dL/dμ, or equivalently, λ = dL/dN, where N is the number of fluid molecules on the pore wall. In the spirit of the mean-field DFT that we utilize, we also estimate ϵff using the mean-field theory that relates the bulk critical temperature Tc to ϵff, namely, ϵff/kB = 4Tc/Z. Thus, with Z = 4, we obtain ϵff/kB = Tc. Therefore, ϵff/kB ≈ 304.25 and 190.82 K for, respectively, CO2 and CH4. Although as discussed earlier, all parameters can be computed by direct MD simulations; in this paper, we estimate the other two parameters, ϵfm and λ, by fitting the model to the available experimental data at various temperatures. As we show below, the two parameters are not very sensitive to temperature. The resulting estimates of the two parameters are presented in Tables 1 and 2. Note that negative and positive values of λ correspond, respectively, to swelling and contraction of the porous materials. We will come back to these points shortly. 4.2. Adsorption Isotherms. The amount adsorbed is calculated by ρ=

1 ∑ρ M i i

T (K)

ϵ* (np)

λ* (np)

ϵ* (lp)

λ* (lp)

254 273 298 343

3 2.9 2.7

20 20 20

2.2 1.9 1.6 1.35

−30 −30 −30 −30

Table 2. Dependence of the Parameters λ* and ϵ* on Temperature for CH4 in the np and lp Phases T (K)

ϵ* (np)

λ* (np)

ϵ* (lp)

λ* (lp)

183 196 213 250

2.5 2.4 2.5

20 20 20

2.2 2.0 1.9 1.8

−30 −30 −30 −30

of lattice units and occupancy. In order to compare the computed isotherms with the experimental data, we must convert the units to physical ones. This is done through the maximum adsorption capacity, nmax, and Vp, with nmax and Vp being, respectively, the maximum density of the adsorbed phase and the sorption pore volume. Theoretically, VvdW−1 ≤ nmax ≤ nϕ, where VvdW is the van der Waals covolume and nϕ is the number density of a close packing of gas particles. For CO2, for example, the lower bound agrees with experimental data for nmax. Thus, because the maximum adsorption capacity in each phase of the flexible MOF that we study has already been reported,43 we utilize it in our computations to convert the computed isotherms to physical units. To compute the adsorption isotherms, the bulk chemical potential μb is increased incrementally from a very large negative value to the saturation bulk chemical potential, −Zϵff/ 2. To convert the chemical potential to pressure, the standard mean-field equations of state of a lattice gas are utilized, which are given by44 1 ϕpN = −NkBT ln(1 − ρb ) − ZNεff ρb 2 (11) 2 ⎛ ρ ⎞ b ⎟⎟ − Zεff ρb μ b = kBT ln⎜⎜ ⎝ 1 − ρb ⎠

(12)

where ϕp = PVs is the lattice equivalent of the bulk pressure, ρb is the density of bulk fluid, and Vs is the physical volume of (the equivalent to) a lattice site. Figure 3 compares the computed adsorption isotherms of CO2 in MIL-53(Al) with the corresponding experimental data at the same temperatures. The isotherms at 254, 273, and 298 K undergo a jump, slightly similar to a first-order (discontinuous) phase transition, which is typical for type IV isotherms45,46 that mesoporous materials exhibit. The sharp transition for intermediate values of pressure emanates from the aforementioned breathing phenomenon and indicates a structural transition between the np and lp phases. When the isotherm is of type IV, the np phase is the most energetically stable form, and to compute the isotherms, the empty framework prior to loading it is in that phase. The computed isotherms for CH4 are presented in Figure 4, where they are compared with the experimental data.43 Similar to CO2, the isotherms at the three lowest temperatures are of type IV, which is characterized by the large jump. On the other hand, the CO2 isotherm at 343 K and that of CH4 at 250 K are

(10)

where M is the number of lattice sites in the pore space. The adsorbed amount or fluid density calculated by eq 10 is in terms 9468


Article


Figure 3. Comparison of the computed CO2 adsorption isotherms with the experimental data.

Figure 4. Comparison of the computed CH4 adsorption isotherms with the experimental data.

and 2 for both np and lp phases, where we present ϵ* = ϵfm/ϵff and λ* = λ/ϵff. Several points are worth pointing out. (i) The strain energy parameter λ* is the same for both CO2 and CH4. This is expected because λ* is a characteristic of the porous material. (ii) λ* is negative for the lp phase but positive for the np phase, indicating, as pointed out earlier, swelling of the material in the former phase but contraction in the latter one. This demonstrates the internal consistency of the model.

of type I, according to the classification of the International Union of Pure and Applied Chemistry,44,46 which develops when adsorption is limited to at most a few molecular layers that, naturally, occur at higher temperatures. In this case, it is the lp phase that constitutes the most stable phase. In the context of the theoretical model presented here, each phase, np or lp, is characterized by its own values of two free parameters of the model. We estimated the two fitting parameters and their dependence on the temperature for both CO2 and CH4; they are listed, respectively, in Tables 1 9469


Article


Figure 5. Pressure dependence of the fractions of the np and lp phases in the transition zone during CO2 adsorption.

(iii) The absolute value of λ* for the lp phase is higher than that for the np phase, implying higher volumetric strain in the lp framework. (iv) For both CO2 and CH4, the parameter ϵ* for the np phase varies only weakly with temperature. In fact, ϵ* = 2.85 ± 0.15 represents an accurate average value, with an estimated variation of only 5.2%. Similarly, for CH4, the estimate ϵ* = 2.45 ± 0.05 is a highly accurate estimate. (v) ϵ* for both CO2 and CH4 in the np phase is larger than that for the lp phase, hence suggesting higher affinity of the np framework for adsorption of both CO2 and CH4. This is in agreement with the available data1,2,43 and indicates once again the internal consistency of our formulation. Note also that ϵ* for CO2 in the np phase is higher than that for CH4 in the same phase, indicating that CO2 has a stronger interaction with the framework. (vi) The parameter ϵ* for the lp phase exhibits a larger range of variations, but even in this case, the variations of ϵ* for both gases are only by a factor of about 1.6. (vii) The difference between the values of ϵ* for CO2 and CH4 in both phases is not large. As Figures 3 and 4 indicate, the model reproduces the experimental data for the np phase up to the pressures at which the transition in the isotherms commences. At higher pressures, the computed isotherms at low temperatures exhibit a plateau that corresponds to the maximum amount of adsorbed CO2 and CH4 or the saturation capacity of the np phase.43 The model also provides accurate isotherms in the lp phase at pressures larger than those in the transition zone. It would perhaps be instructive to compare the present model with that of Neimark et al.,37 which is considered to be an accurate model for sorption in and deformation of MOFs. A fitting parameter in the study by Neimark et al.37 is the Henry constant KH, which represents the adsorption affinity of adsorbates for the MOF structure. They found that KH for the np phase is larger than that of the lp phase. In our model, the parameter ϵ* has the same characteristics: it is representative of the adsorption affinity of the gases, and for each adsorbate (CO2 or CH4) at any temperature, ϵ* for the np phase is larger than that for the lp phase. Thus, at a qualitative level, the two parameters are equivalent and could perhaps be related as well. The model proposed by Neimark et al. contains two other parameters, namely, dN0/dVc and dKH/dVc, where N0 is the unit cell capacity and Vc is the volume of the cell. Both parameters influence the adsorption stress. Neimark et al.37 found that, at low pressures, the effect of dKH/dVc (which is

negative) on the stress is more prominent, giving rise to negative of adsorption stress and hence resulting in the contraction of MOFs. At high pressures, on the other hand, positive adsorption stress is exerted on the sample, causing material expansion. In the model proposed here, the parameter λ has essentially the same role as the two parameters of Neimark et al.’s. As described earlier, λ represents the coupling between the fluid and solid phases and is the amount of deformation of the material as a result of adsorption of one unit of adsorbate. In addition, λ is positive in the np phase, indicating a negative adsorption stress and sample contraction, and negative, which is representative of positive adsorption stress and sample expansion, in the lp phase. 4.3. Plastic Deformation and Phase Fractions. Neimark et al.22,37 suggested that the step-wise transition in the isotherms arises from a large plastic deformation in the structure of MOFs. If so, to model the transition zone, one needs to substitute eq 3 for the linearly elastic strain energy with one for the energy stored in a plastic material.47−50 Such an expression depends, of course, on the type of the material under deformation. The rest of the model will remain the same. Work in this direction is in progress, and the results will be reported in a future paper. For now, we attempt to develop an understanding of the root cause of the jumps in the adsorption isotherms that, as mentioned earlier, are representative of a structural phase transition. To do so, we combine our model with a phase-mixture model.50 Previous studies4,22 indicated that there is a mixture of the np and lp phases in the transition zone. Thus, suppose that Xnp and Xlp are the fractions of the two phases in the pressure range of the transition zone. In other words, the hypothesis is that the system in the transition zone is neither pure np nor lp phase but a mixture of the two. Then, the total amount of adsorbed CO2, ρt, is computed by model ρt = X npρnp + Xlpρlpmodel

(13)

where ρmodel and ρmodel are the computed adsorbed amounts for np lp the pure np and lp phases in the same pressure ranges, respectively. Because Xnp + Xlp = 1, eq 13 is rewritten as model ρt = X npρnp + ρlpmodel (1 − X np)

(14)

There are two unknowns, but only one equation, eq 14. Thus, to determine Xnp, the adsorption amount in the pure np and lp phases, that is, ρmodel and ρmodel , are computed in the np lp pressure range that corresponds to the transition zone. Then, Xnp is selected in such a way that ρt matches the experimental data. 9470


Article


Figure 6. Pressure dependence of the fraction of the np and lp phases in the transition zone during adsorption of CH4.

Figure 7. Computed volumetric strain in MIL-53(Al) during CO2 adsorption in the lp and np frameworks.

5. SUMMARY We described a model for sorption-induced deformation of MOF materials. Although the model was used to simulate and compute adsorption of CO2 and CH4 in flexible MIL-53(Al) MOFs and their associated induced deformation, it is, as discussed earlier, completely general and may be applicable to a wide variety of other types of materials, including other subclasses of MOF materials in which deformation is induced by the interactions between the material and the guest molecules. The important point to emphasize is that, although we assumed that the MOF materials follow linear elastic deformation, other types of deformation, including plastic deformation in the transition zone, may also be considered within the framework of the model, as long as the constitutive equation that describes the deformation is available. The model produces all important aspects of sorption of CO2 and CH4 in MIL-53(Al), including distinct features of the well-known np and lp phases, coexistence of the two phases at intermediate pressures that is related to a structural transition, and contraction and expansion of the material at low and high pressures and temperatures.

Figures 5 and 6 present, respectively, the computed Xnp and Xlp as a function of the pressure in the transition zone for CO2 and CH4. They indicate that the fractions of the np and lp phases in the transition region in which the breathing phenomenon occurs are strong functions of the pressure. Thus, the transition is indeed a manifestation of a structural transition in the MOF. 4.4. Strain and Deformation. As described earlier, numerical simulation of the model yields both the sorption isotherms and the displacement field throughout the material, resulting from the elastic deformation. Figures 7 and 8 present, respectively, the volumetric strain “isotherms” for the np and lp phases upon adsorption of CO2 and CH4. During adsorption of the gases, the strain in the np phase at the three lowest temperatures is negative, which is indicative of the framework contraction. At the highest temperature, 343 K for CO2 and 250 K for CH4, the strain is positive and corresponds to the expansion of the framework in the lp phase. These results are consistent not only with those shown in Figures 3 and 4 but also with the available experimental data.43 9471


Article


Figure 8. Computed volumetric strain in MIL-53(Al) during CH4 adsorption in the lp and np frameworks.

■

(7) Perry, J. J.; Perman, J. A.; Zaworotko, M. J. Design and Synthesis of Metal-Organic Frameworks Using Metal-Organic Polyhedra as Supermolecular Building Blocks. Chem. Soc. Rev. 2009, 38, 1400− 1417. (8) Wilmer, C. E.; Leaf, M.; Lee, C. Y.; Farha, O. K.; Hauser, B. G.; Hupp, J. T.; Snurr, R. Q. Large-Scale Screening of Hypothetical MetalOrganic Frameworks. Nat. Chem. 2011, 4, 83−89. (9) Chung, Y. G.; Camp, J.; Haranczyk, M.; Sikora, B. J.; Bury, W.; Krungleviciute, V.; Yildirim, T.; Farha, O. K.; Sholl, D. S.; Snurr, R. Q. Computation-Ready, Experimental Metal-Organic Frameworks: A Tool to Enable High-Throughput Screening of Nanoporous Crystals. Chem. Mater. 2014, 26, 6185−6192. (10) Bae, Y.-S.; Yazaydin, A. Ö .; Snurr, R. Q. Evaluation of the BET Method for Determining Surface Areas of MOFs and Zeolites that Contain Ultra-Micropores. Langmuir 2010, 26, 5475−5483. (11) Gómez-Gualdrón, D. A.; Colón, Y. J.; Zhang, X.; Wang, T. C.; Chen, Y.-S.; Hupp, J. T.; Yildirim, Y.; Farha, O. K.; Zhang, J.; Snurr, R. Q. Evaluating Topologically Diverse Metal-Organic Frameworks for Cryo-Adsorbed Hydrogen Storage. Energy Environ. Sci. 2016, 9, 3279− 3289. (12) Li, H.; Eddaoudi, M.; O’Keeffe, M.; Yaghi, O. M. Design and Synthesis of an Exceptionally Stable and Highly Porous Metal-Organic Framework. Nature 1999, 402, 276−279. (13) Chen, B.; Eddaoudi, M.; Hyde, S. T.; O’Keeffe, M.; Yaghi, O. M. Interwoven Metal-Organic Framework on a Periodic Minimal Surface with Extra-Large Pores. Science 2001, 291, 1021−1023. (14) Furukawa, H.; Ko, N.; Go, Y. B.; Aratani, N.; Choi, S. B.; Choi, E.; Yazaydin, A. O.; Snurr, R. Q.; O’Keeffe, M.; Kim, J.; et al. UltraHigh Porosity in Metal-Organic Frameworks. Science 2010, 329, 424− 428. (15) Mishra, P.; Edubilli, S.; Mandal, B.; Gumma, S. Adsorption of CO2, CO, CH4 and N2 on DABCO Based Metal-Organic Frameworks. Microporous Mesoporous Mater. 2013, 169, 75−80. (16) Mason, J. A.; Sumida, K.; Herm, Z. R.; Krishna, R.; Long, J. R. Evaluating Metal-Organic Frameworks for Post-Combustion Carbon Dioxide Capture via Temperature Swing Adsorption. Energy Environ. Sci. 2011, 4, 3030−3040. (17) Finsy, V.; Ma, L.; Alaerts, L.; De Vos, D. E.; Baron, G. V.; Denayer, J. F. M. Separation of CO2/CH4 Mixtures with the MIL53(Al) Metal-Organic Framework. Microporous Mesoporous Mater. 2009, 120, 221−227.

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Phone: (213) 740-2064. ORCID

Muhammad Sahimi: 0000-0002-8009-542X Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS S.B. is grateful to the Mork Family Department of Chemical Engineering and Materials Science at the University of Southern California for a two-year doctoral fellowship. This research was supported as part of the Center for Geologic Storage of CO2, an Energy Frontier Research Center funded by the U.S. Department of Energy (DOE), Office of Science, Basic Energy Sciences (BES), under Award number DE-SC0C12504. The calculations were carried out using the computer network of the High Performance Computing Center at the University of Southern California.

■

REFERENCES

(1) Alhamami, M.; Doan, H.; Cheng, C.-H. A Review on Breathing Behaviors of Metal-Organic Frameworks (MOFs) for Gas Adsorption. Materials 2014, 7, 3198−3250. (2) Mishra, P.; Uppara, H. P.; Mandal, B.; Gumma, S. Adsorption and Separation of Carbon Dioxide Using MIL-53(Al) Metal-Organic Framework. Ind. Eng. Chem. Res. 2014, 53, 19747−19753. (3) Special issue of Chem. Soc. Rev. on MOFs, 2014, 43. (4) Ghysels, A.; Vanduyfhuys, L.; Vandichel, M.; Waroquier, M.; van Speybroeck, V.; Smit, B. On the Thermodynamics of Framework Breathing: A Free Energy Model for Gas Adsorption in MIL-53. J. Phys. Chem. C 2013, 117, 11540−11554. (5) Yaghi, O. M.; O’Keeffe, M.; Ockwig, N. W.; Chae, H. K.; Eddaoudi, M.; Kim, J. Reticular Synthesis and the Design of New Materials. Nature 2003, 423, 705−714. (6) Férey, G. Hybrid Porous Solids: Past, Present, Future. Chem. Soc. Rev. 2008, 37, 191−214. 9472


Article

The Journal of Physical Chemistry C (18) Alaerts, L.; Maes, M.; Giebeler, L.; Jacobs, P. A.; Martens, J. A.; Denayer, J. F. M.; Kirschhock, C. E. A.; De Vos, D. E. Selective Adsorption and Separation of ortho-Substituted Alkylaromatics with the Microporous Aluminum Terephthalate MIL-53. J. Am. Chem. Soc. 2008, 130, 14170−14178. (19) Herm, Z. R.; Swisher, J. A.; Smit, B.; Krishna, R.; Long, J. R. Metal-Organic Frameworks as Adsorbents for Hydrogen Purification and Precombustion Carbon Dioxide Capture. J. Am. Chem. Soc. 2011, 133, 5664−5667. (20) Beurroies, I.; Boulhout, M.; Llewellyn, P. L.; Kuchta, B.; Férey, G.; Serre, C.; Denoyel, R. Using Pressure to Provoke the Structural Transition of Metal-Organic Frameworks. Angew. Chem., Int. Ed. 2010, 49, 7526−7529. (21) Liu, Y.; Her, J.-H.; Dailly, A.; Ramirez-Cuesta, A. J.; Neumann, D. A.; Brown, C. M. Reversible Structural Transition in MIL-53 with Large Temperature Hysteresis. J. Am. Chem. Soc. 2008, 130, 11813− 11818. (22) Neimark, A. V.; Coudert, F.-X.; Triguero, C.; Boutin, A.; Fuchs, A. H.; Beurroies, I.; Denoyel, R. Structural Transitions in MIL-53 (Cr): View from Outside and Inside. Langmuir 2011, 27, 4734−4741. (23) Bourrelly, S.; Llewellyn, P. L.; Serre, C.; Millange, F.; Loiseau, T.; Férey, G. Different Adsorption Behaviors of Methane and Carbon Dioxide in the Isotypic Nanoporous Metal Terephthalates MIL-53 and MIL-47. J. Am. Chem. Soc. 2005, 127, 13519−13521. (24) Kitagawa, S.; Uemura, K. Dynamic Porous Properties of Coordination Polymers Inspired by Hydrogen Bonds. Chem. Soc. Rev. 2005, 34, 109−119. (25) Serre, C.; Mellot-Draznieks, C.; Surble, S.; Audebrand, N.; Filinchuk, Y.; Ferey, G. Role of Solvent-Host Interactions That Lead to Very Large Swelling of Hybrid Frameworks. Science 2007, 315, 1828−1831. (26) Mellot-Draznieks, C.; Serre, C.; Surblé, S.; Audebrand, N.; Férey, G. Very Large Swelling in Hybrid Frameworks: A Combined Computational and Powder Diffraction Study. J. Am. Chem. Soc. 2005, 127, 16273−16278. (27) Mota, J. P. B.; Martins, D.; Lopes, D.; Catarino, I.; Bonfait, G. Structural Transitions in the MIL-53(Al) Metal-Organic Framework upon Cryogenic Hydrogen Adsorption. J. Phys. Chem. C 2017, 121, 24252−24263. (28) Moghadam, P. Z.; Ivy, J. F.; Arvapally, R. K.; dos Santos, A. M.; Pearson, J. C.; Zhang, L.; Tylianakis, E.; Ghosh, P.; Oswald, I. W. H.; Kaipa, U.; et al. Adsorption and Molecular Siting of CO2, Water, and other Gases in the Superhydrophobic, Flexible Pores of FMOF-1 from Experiment and Simulation. Chem. Sci. 2017, 8, 3989−4000. (29) Bon, V.; Klein, N.; Senkovska, I.; Heerwig, A.; Getzschmann, J.; Wallacher, D.; Zizak, I.; Brzhezinskaya, M.; Mueller, U.; Kaskel, S. Exceptional Adsorption-Induced Cluster and Network Deformation in Flexible Metal-Organic Framework DUT-8 (Ni) Observed by In-Situ X-Ray Diffraction and EXAFS. Phys. Chem. Chem. Phys. 2015, 17, 17471−17479. (30) Simon, C. M.; Braun, E.; Carraro, C.; Smit, B. Statistical Mechanical Model of Gas Adsorption in Porous Crystals with Dynamic Moieties. Proc. Natl. Acad. Sci. U.S.A. 2017, 114, E287−E296. (31) Reichenbach, C.; Kalies, G.; Lincke, J.; Lässig, D.; Krautscheid, H.; Moellmer, J.; Thommes, M. Unusual Adsorption Behavior of a Highly Flexible Copper-Based MOF. Microporous Mesoporous Mater. 2011, 42, 592−600. (32) Tanaka, H.; Hiraide, S.; Kondo, A.; Miyahara, M. T. Modeling and Visualization of CO2 Adsorption on Elastic Layer-Structured Metal-Organic Framework-11: Toward a Better Understanding of Gate Adsorption Behavior. J. Phys. Chem. C 2015, 119, 1533−1543. (33) Wu, H.; Thibault, C. G.; Wang, H.; Cychosz, K. A.; Thommes, M.; Li, J. Effect of Temperature on Hydrogen and Carbon Dioxide Adsorption Hysteresis in an Ultramicroporous MOF. Microporous Mesoporous Mater. 2016, 219, 186−189. (34) Yang, Q.; Liu, D.; Zhong, C.; Li, J.-R. Development of Computational Methodologies for Metal-Organic Frameworks and Their Application in Gas Separations. Chem. Rev. 2013, 113, 8261− 8323.

(35) Yue, Y.; Rabone, J. A.; Liu, H.; Mahurin, S. M.; Li, M.-R.; Wang, H.; Lu, Z.; Chen, B.; Wang, J.; Fang, Y.; Dai, S. A Flexible MetalOrganic Framework: Guest Molecules Controlled Dynamic Gas Adsorption. J. Phys. Chem. C 2015, 119, 9442−9449. (36) Coudert, F.-X.; Jeffroy, M.; Fuchs, A. H.; Boutin, A.; MellotDraznieks, C. Thermodynamics of Guest-Induced Structural Transitions in Hybrid Organic-Inorganic Frameworks. J. Am. Chem. Soc. 2008, 130, 14294−14302. (37) Neimark, A. V.; Coudert, F.-X.; Boutin, A.; Fuchs, A. H. StressBased Model for the Breathing of Metal-Organic Frameworks. J. Phys. Chem. Lett. 2010, 1, 445−449. (38) Pitard, E.; Rosinberg, M. L.; Stell, G.; Tarjus, G. Critical Behavior of a Fluid in a Disordered Porous Matrix: an OrnsteinZernike Approach. Phys. Rev. Lett. 1995, 74, 4361−4364. (39) Kierlik, E.; Rosinberg, M. L.; Tarjus, G.; Pitard, E. MeanSpherical Approximation for a Lattice Model of a Fluid in a Disordered Matrix. Mol. Phys. 1998, 95, 341−351. (40) Kierlik, E.; Monson, P. A.; Rosinberg, M. L.; Sarkisov, L.; Tarjus, G. Capillary Condensation in Disordered Porous Materials: Hysteresis versus Equilibrium Behavior. Phys. Rev. Lett. 2001, 87, 055701. (41) Guyer, R. A.; Kim, H. A. Theoretical Model for Fluid-Solid Coupling in Porous Materials. Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys. 2015, 91, 042406. (42) Bakhshian, S.; Sahimi, M. Adsorption-Induced Swelling of Porous Media. Int. J. Greenhouse Gas Control 2017, 57, 1−13. (43) Boutin, A.; Coudert, F.-X.; Springuel-Huet, M.-A.; Neimark, A. V.; Férey, G.; Fuchs, A. H. The Behavior of Flexible MIL-53(Al) upon CH4 and CO2 Adsorption. J. Phys. Chem. C 2010, 114, 22237−22244. (44) Monson, P. A. Understanding Adsorption/Desorption Hysteresis for Fluids in Mesoporous Materials Using Simple Molecular Models and Classical Density Functional Theory. Microporous Mesoporous Mater. 2012, 160, 47−66. (45) Sing, K. S. W. Reporting Physisorption Data for Gas/Solid Systems with Special Reference to the Determination of Surface Area and Porosity. Pure Appl. Chem. 1985, 57, 603−619. (46) Sahimi, M. Flow and Transport in Porous Media and Fractured Rock, 2nd ed.; Springer: New York, 2011; Chapter 5. (47) Lubliner, J. Plasticity Theory; Macmillan: New York, 1990. (48) Hill, R. The Mathematical Theory of Plasticity; Oxford University Press: London, 1998. (49) Computational Plasticity; Onate, E., Owen, R., Eds.; Springer: Dordrecht, 2007. (50) Ghoufi, A.; Maurin, G. Hybrid Monte Carlo Simulations Combined with a Phase Mixture Model to Predict the Structural Transitions of a Porous Metal-Organic Framework Material upon Adsorption of Guest Molecules. J. Phys. Chem. C 2010, 114, 6496− 6502.

9473


Theoretical Model and Numerical Simulation of Adsorption and

Recommend Documents