Thermodynamic Route to Efficient Prediction of Gas Diffusivity in

Sep 15, 2017 - accurate prediction of the self-diffusivity of gas molecules in nanoporous materials ... classical molecular dynamics (MD) simulation.5...
6 downloads 0 Views 1MB Size
Subscriber access provided by UNIV OF CAMBRIDGE

Article

A Thermodynamic Route to Efficient Prediction of Gas Diffusivity in Nanoporous Materials Yun Tian, Xiaofei Xu, and Jianzhong Wu Langmuir, Just Accepted Manuscript • DOI: 10.1021/acs.langmuir.7b02428 • Publication Date (Web): 15 Sep 2017 Downloaded from http://pubs.acs.org on September 20, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Langmuir is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 24

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

A Thermodynamic Route to Efficient Prediction of Gas Diffusivity in Nanoporous Materials Yun Tian1, Xiaofei Xu2 and Jianzhong Wu1* 1

Department of Chemical and Environmental Engineering, University of California, Riverside, CA 92521, United States

2

Center for Soft Condensed Matter Physics and Interdisciplinary Research, Soochow University, Suzhou 215006, China Abstract

We report an efficient computational procedure for rapid and accurate prediction of the self-diffusivity of gas molecules in nanoporous materials by implementing the transition state theory for inter-cage hopping at infinite dilution with the string method in conjunction with the excess-entropy scaling for predicting gas diffusion coefficients at finite loadings. The theoretical procedure has been calibrated with molecular dynamics simulations for the diffusion coefficients of methane and hydrogen gases in representative nanoporous materials including metal organic frameworks and zeolites. Combined with the classical density functional theory for calculating the excess entropy of gas molecules in micropores, the theoretical procedure enables efficient computation of both thermodynamic and transport properties important for design and screening of nanostructured materials for gas storage and separation.

*

To whom correspondence should be addressed. Email: [email protected] 1 ACS Paragon Plus Environment

Langmuir

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 24

1. Introduction Molecular transport in nanoporous materials plays a crucial role in emerging nanotechnologies for gas storage and separation1. Although practical applications are often concerned

with

transport

or

Fickian

diffusion

coefficients

affiliated

with

various

phenomenological equations, self-diffusivity of gas molecules in nanoporous materials has been widely studied through both experimental and theoretical means2, 3. The quantity is significant not only from a practical perspective but also for theoretical developments. Self-diffusivity is closely connected with different types of transport coefficients and can be directly measured from pulsed-field gradient-nuclear magnetic resonance (PFG-NMR) experiment4. Theoretical predictions of self-diffusivity are mostly based on classical molecular dynamics (MD) simulation5. For gas molecules in nanoporous materials, the self-diffusion coefficient can be obtained by analyzing the trajectories of gas molecules according to Einstein’s equation. While the computational procedure is straightforward and generically applicable to all nanoporous systems, one of the major difficulties with direct MD simulation is that the selfdiffusion coefficients of gas molecules in a nanostructured material are typically very small (~ 10 −12 m 2 / s ).

The slow diffusion process is often beyond the time scale accessible to

conventional simulation methods6. Among various indirect procedures for calculating selfdiffusivity, the transition state theory (TST) represents a popular choice. TST predicts kinetic coefficients of infrequent events such as gas hopping among various cages of nanoporous materials7-10. The thermodynamic approach can be orders of magnitude more efficient than “brute force” MD simulations while still retaining full atomistic details. Numerical implementation of TST for predicting gas diffusivity requires efficient algorithms to evaluate the diffusion pathways and the underlying free-energy landscape (viz.,

2 ACS Paragon Plus Environment

Page 3 of 24

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

reaction coordinate diagram). The free-energy landscape defines diffusion paths for a gas molecule hopping from the original state to the destination state. For gas diffusion at infinite dilution within a rigid adsorbent structure, the free energy is replaced by a one-body potential that arises from the interaction of the gas molecule with the adsorbent atoms. In that case, one may implement TST by an exhaustive enumeration of the energy landscape over the entire volume11. While the numerical procedure is in principle exact, construction of the reaction coordinate diagram is computationally demanding and not easily extendable to polyatomic molecules and to free-energy calculations as required for predicting diffusivity at finite loadings12. In the latter case, TST is often implemented by using a predefined reaction coordinate that is obtained from the pore geometry of a specific porous material under consideration3. In combination with the Bennett–Chandler theory for dynamic corrections to account for the probability of recrossing at the transition state, the simplified method was found in good agreement with direct MD simulation at arbitrary loadings13. To a certain degree, the dynamic correction, which also relies on MD simulation, compensates numerical errors introduced in locating the transition state and the reaction coordinate.14,

15

Regrettably, computation of the

dynamic coefficients is itself nontrivial, in particular when the reaction pathway involves multiple pores16. Given a potential (or free) energy landscape, the string method is emerged as one of the most efficient ways to estimate the multidimensional reaction coordinate of diverse kinetic processes17. In contrast to conventional methods that rely on grid searching and cluster analysis11, the string method is able to identify the diffusion pathways based on the energy gradients such that each path follows an exactly minimum energy route. Besides, the string method does not require a predefined reaction coordinate and is directly applicable to systems with complicated

3 ACS Paragon Plus Environment

Langmuir

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

curvilinear paths18. Whereas the string method and its later modifications have been extensively used in studying the kinetics of various rare events, we are unaware of previous publication on its application to gas diffusion in nanoporous materials. The scaling relation between self-diffusivity and excess entropy was discovered first by Rosenfeld on the basis of earlier simulation results for the transport coefficients of a wide variety of bulk systems19. In our previous work20, we extend Rosenfeld’s excess-entropy scaling method to gas diffusion in nanoporous materials by combination with the Knudsen model. The basic idea is to divide self-diffusivity into contributions due to adsorbate-adsorbate and adsorbent-adsorbate interactions. In this work, we implement TST in conjunction with the string method to define the reaction coordinate of gas diffusion in nanoporous materials at infinite dilution. As reported in our earlier publication20, the diffusion coefficients at finite loadings can be predicted from excess-entropy scaling in combination with the classical density functional theory (cDFT) calculations. We test the new procedure with MD simulation for the diffusion of light gas molecules (H2 and CH4) in several nanoporous materials including metal organic frameworks (MOF) and zeolites. By eschewing a priori definition of reaction coordinates and direct evaluation of the potential energy landscape, the theoretical procedure is found more robust than existing methods for predicting gas diffusivity at infinite dilution. 2. Theoretical Methods String method To illustrate the basic concepts, consider first the diffusion of a light gas molecule in a nanoporous material at infinite dilution. We assume that the gas molecule has a spherical shape and the adsorbent has a rigid structure. The potential energy landscape is defined by the external energy of the gas molecule, V ext (r) , due to its interaction with the adsorbent atoms

4 ACS Paragon Plus Environment

Page 4 of 24

Page 5 of 24

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

V ext (r) = ∑ ugi (r − ri )

(1)

i

where r stands for the position of the gas molecule (subscript g), and ugi is its pair interaction energy with atom i of the adsorbent located at ri . For simplicity, we assume that the pair potential is given by the Lennard-Jones (LJ) model

 σ  12  σ  6  ugi (r) = 4ε gi  gi  −  gi    r   r    

(2)

where ε gi and σ gi are the energy and size parameters, respectively. For all materials considered in this work, the adsorbent structures are obtained from experimental Crystal Information Files (CIFs)21,

22

and the LJ parameters are from the UFF force field for each atom with the

conventional combining rules for different atoms23. The same force field parameters are used in MD simulations to benchmark the theoretical development. The one-body potential allows us to identify the positions of minimum energy corresponding to the initial state and the final state of gas hopping between neighboring cages of the nanoporous material. In contrast to conventional methods that define the reaction coordinate in terms of a straight line along the pore centers, we identify a multidimensional minimumenergy path by implementing the string method. Specifically, we define s as a normalized reaction coordinate, i.e., 0 ≤ s ≤ 1 , and p(s) as a minimum energy pathway that connects the initial ( s = 0 ) and final ( s = 1 ) states of gas molecule hopping. The saddle point along the diffusion pathway, which has a maximum energy and is close to the void center linking neighboring cages, is identified as the transition state. The diffusion pathway is calculated in such a way that tangent vector

along the

minimum energy path (dimensionless) is parallel to the direction of energy gradient 5 ACS Paragon Plus Environment

Langmuir

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 24

∇V ≡ ∇V ext / | ∇V ext | . In other words, the component of the energy gradient normal to the diffusion pathway vanishes. Mathematically, the tangent vector and the (dimensionless) energy gradient satisfies17 , (3) where

(4) is the normalized tangent along the diffusion pathway, and * donates the inner product defined as f * g = ∫ f ( r ) g ( r ) dr . The diffusion pathway, p(s) , shows a minimum energy along with

reaction coordinate s(r) and can be identified from Eq.(3) and (4) following a modified steepestdescent algorithm24, 25 (5) where Λ denotes the Lagrange multiplier, and t is a fictitious time (it has the units of length) for evolving the equations on the potential energy landscape. The Lagrange multiplier is introduced to normalize the arc length. For numerical integrations, the grid size is fixed at 0.1 Å for all cases considered in this work. It is worth noting that the string method allows us to identify a reaction coordinate s (r ) that accounts for the 3-dimensional structure of nanopores including those with straight and zigzag geometries. The procedure is directly applicable to systems with other forms of gasadsorbent interactions (e.g., Coulomb plus LJ interactions). By replacing the one-body potential energy with the local excess chemical potential inside the nanoporous material, the string method can be similarly applied to non-spherical gas molecules even when the adsorbent structure is not 6 ACS Paragon Plus Environment

Page 7 of 24

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

fixed. The main difference in implementation of the string method is that, in general, the excess chemical potential calculation requires statistical-mechanical calculations such as simulation methods12 or classical density functional theory (cDFT)26. Transition state theory We use the transition state theory (TST) to describe the kinetics of gas hopping from one cage to another in a nanoporous material. At infinite dilution, the self-diffusion coefficient is related to the transmission by

D0 =

1 2 ka 2

(6)

where a is the distance between the gas molecule in the neighboring cages (the initial and final states of the transmission), and k is the hopping rate. According to a simplified version of TST10, the hoping rate is determined by

k = v exp(−∆ E/ k BT )

(7)

where ∆E is the potential energy barrier in the diffusion path, k B is the Boltzmann constant, T is the absolute temperature, and v is the transmission frequency. The potential energy barrier is defined as the difference between the gas molecule in the transition state and that in the initial state. In this work, the pre-exponential factor is taken as v = 1012 s-1 for metal-organic frameworks (MOFs)10 and v = 1011 s-1 for zeolites27. These values are not sensitive to the chemical type of the adsorbents within the same category of nanoporous materials. The simplified TST method is valid only when the potential energy barrier is sufficiently high compared with other barriers along the diffusion pathway. For materials with arbitrary pore structures, the hopping rate is calculated from the following more general equation28

7 ACS Paragon Plus Environment

Langmuir

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

k BT k= 2π m

exp  − V ext (s* ) / k B T 



1

0

exp  − V ext (s) / k B T  ds

Page 8 of 24

(8)

where m stands for the mass of the gas molecule, the integral is performed along the reaction coordinate, and superscript * represents the transition state. Because this work is only concerned with diffusion of light-gas molecules such as hydrogen and methane in rigid nanoporous materials, we expect that the dynamic correction is relatively insignificant10, 27. Like the string method, TST is applicable to gas diffusion at both infinite dilution and finite loadings. In principle, it is able to account for the flexibility of the adsorbent structure as well as the dynamic nature of gas molecules14, 29. Whereas flexibility may play a major role on dynamic processes in general18, it has been shown that its influence is relatively insignificant on the diffusivity of light gas molecules in a rigid nanostructured material30. Excess entropy scaling One major advantage in application of TST at infinite dilution is that the energy landscape can be directly calculated from the one-body external potential. At high loadings, we need to account for gas-gas interactions even for rigid framework materials and computation of the local free energy is in general computationally demanding. In our previous publications20, 31, we demonstrate that the self-diffusivity of a gas molecule at finite loadings can be predicted from a linear combination of that at the infinite dilution and excess due to interaction between gas molecules:

 ασ 3  ασ 3 ln Ds =  1− ln D + ln DE  0 v free  v free 

(9)

where D0 is the self-diffusivity at the infinite dilution, α is a constant related to the maximum molecular packing density, v free represents the total accessible volume per gas molecule in the 8 ACS Paragon Plus Environment

Page 9 of 24

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

porous material. The mixing parameter, ασ 3 / v free accounts for the relative contribution of gasgas interactions to the self-diffusivity of gas molecules. At maximum molecular packing density, the self-diffusivity is fully determined by the interactions between gas molecules. Therefore, parameter α may be estimated from the close packing limit of spherical particles, α = 2 / 2 . While this parameter changes slightly for different systems depending on the crystalline structure, it is independent of the thermodynamic conditions of the bulk gas including pressure and temperature20. As discussed above, we calculate D0 from TST implemented with the string method. The contribution due to gas-gas interactions, DE , is calculated from the excess entropy scaling method20:

DE =

(

0.585 kBT exp 0.788sex / k B 1/3 m ρ

)

(10)

where ρ is the average number density of gas molecules inside the porous material, s ex is the excess entropy per confined gas molecule calculated from the classical density functional theory (cDFT)26. 3. Results and Discussion The numerical performance of the transition state theory (TST) is dependent on the reaction coordinate diagram, i.e., the variation of the free energy of the entire system during the kinetic process. In this work, we apply TST only to diffusion of a spherical gas molecule at infinite dilution. In that case, the free energy corresponds to a one-body external potential due to interaction of the gas molecule with a rigid porous material. As mentioned above, the procedure can be generalized by using a simulation method to calculate the free energy of a gas molecule in a flexible nanostructure material. 9 ACS Paragon Plus Environment

Langmuir

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 1a) presents the reaction coordinate diagram for CH4 diffusion along the major pore of MOF-5 calculated from the string method. For comparison, we show a similar diagram obtained from the conventional geometric analysis, i.e., assuming that the reaction coordinate coincides with the axis along the pore center. The differences in the free-energy profiles are rather astonishing. First of all, the energy along the diffusion pathway is much smaller than that in the pore center, indicating that the geometric analysis does not provide a faithful description of the minimum energy diffusion pathway. While the external potential shows a smooth sinusoidal curve along the pore center of MOF-5 in harmony with the crystalline structure, its variation along the reaction coordinate is not nearly as smooth, showing some wiggling changes due to the interaction of the methane molecule with individual atoms of the framework. Importantly, the string method allows us to identify multiple metastable states along the reaction pathway for gas diffusion in typical nanoporous materials. The drastic deviation of the diffusion pathway from the pore center contrasts an earlier study for hydrogen diffusion in similar MOF materials32. The string method identifies multiple minimum and maximum energy states along the diffusion pathway, implying that the hopping rate cannot be simply determined from a single rate-limiting barrier as shown in Eq.(7). Instead, we need to calculate the overall diffusion barrier using the more generic equation, Eq.(8). Figure 1b) shows the methane locations corresponding to the major minimum/maximum energy states as identified from Figure 1a). States 1 and 7 identified by the string method are positioned at the most favorable binding (lowest energy) sites of MOF-5 for the methane molecule. These two states define the initial and final states in the transition state theory. The minimum energy states identified from the string methods are distinctively different from those predicted by the geometric analysis. While the latter are located at the pore center, the string method predicts that the minimum energy states are close to the

10 ACS Paragon Plus Environment

Page 10 of 24

Page 11 of 24

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

metal ions (Zn) at the edge of the framework. States 3 and 5 identified from the string method represent the second strongest binding sites of the material, also located close to the metal ions but after the methane molecule passing through two narrow pores. The molecular energy at the center of these narrow pores, which are designated as states 2 and 6 according to the string method, are nearly as high as that corresponding to the transition state (state 4). Figure 1b) shows that the transition state is geometrically accessible; it has a maximum energy along the diffusion pathway because its interaction with the framework material is weakest. In other words, the potential energy barrier for the diffusion of a light gas molecule in nanoporous materials may not arise from the excluded volume effects but from its strong attraction with the binding sites. Because the lower energy states are not located at the pore center as assumed in geometric analysis but near strong binding sites, a priori assignment of the reaction coordinate according to the pore geometry may introduce significant errors to TST calculations. We use the string method to construct the reaction coordinate diagrams for methane and hydrogen diffusion in six representative MOF materials. Figure 2 summarizes the self-diffusion coefficients from the TST predictions in comparison with the corresponding results from direct MD simulation20, 32. Based on the string method, Eq.(6) and Eq.(8), we are able to predict the self-diffusivity of small gas molecules (either CH4 or H2) in excellent agreement with MD data. By contrast, previous application of TST compares poorly with the simulation results as shown by the red circles10,

32

. Without a reliable reaction coordinate diagram, TST can either

overestimate or underestimate the self-diffusivity at infinite dilution. We have also applied the string method for methane diffusion in five zeolites. Figure 3 shows a comparison of the theoretical results from this work with MD data reported by Haldoupis and coworkers27. Again, the numerical performance is much improved with the string

11 ACS Paragon Plus Environment

Langmuir

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 24

method in comparison to previous TST predictions27. For CFI and MTW zeolites, we see noticeable difference between the TST and MD results. The discrepancy is probably due to the uncertainties in MD simulations because sufficient sampling is difficult for calculating diffusivity at the infinite dilution. As demonstrated in our previous publication20, we can predict the gas self-diffusivity at finite loadings using the diffusion coefficient at infinite dilution based on the excess entropy scaling method. Toward that end, we employ the classical density functional theory (cDFT) to predict the adsorption isotherms and the thermodynamic properties of confined gas molecules. As reported in our earlier work20, 31, we formulate the cDFT for inhomogeneous LJ fluids using the weighted density approximation for the Helmholtz energy functional due to the van der Waals attractions and the modified fundamental measure theory to account for the molecular volume exclusive effects. Figure 4 compares the theory and simulation results for the selfdiffusion coefficients of methane in MOF-5 and in CuBTC over a wide range of pressures. We see that the theoretical results are almost within the error bars of direct MD simulation. The good agreement between excess entropy scaling and MD simulation are consistent with the previous calculations for finite gas loadings20. While Figure 4 shows results for gas pressure only up to 100 bars, we expect that the theoretical performance will be satisfactory at higher pressure because, as shown in Eq.(9), at high loading self-diffusivity is primarily determined by the excess-entropy scaling. Indeed, our earlier work shows that the excess entropy scaling method is reliable up to 550 bars. It is worth noting that that the computational time for TST implemented with the string method is generally within 1/3 CPU hour per point (3.0 GHz), while the MD simulation requires at least 10 CPU hours for a single point. Although we are able to implement TST at finite loadings based on the excess chemical potential of gas molecules from cDFT

12 ACS Paragon Plus Environment

Page 13 of 24

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

calculations, the excess entropy scaling is even more convenient because it takes virtually no extra computational cost beyond the adsorption isotherm calculations. Finally, it should be indicated that cDFT calculations provide not only the excess entropy but also other thermodynamic properties including heat effects and adsorption isotherms. For example, Figure 5 shows the adsorption isotherms for methane in MOF-5 and CuBTC at 298 K. As reported in our earlier work31, cDFT predictions agree well with the results from standard GCMC simulation. The good agreement justifies in part the applicability of the excess entropy scaling method to predict the gas diffusivity in nanoporous materials at finite loadings. 4. Conclusion Numerical implementation of the transition state theory (TST) requires a priori knowledge of a free-energy (or potential energy) landscape that defines the gas molecule in the initial and final states as well as a hypersurface that divides the configurational space into the two thermodynamic states. In this work, we demonstrate that the string method provides a more faithful description of the diffusion pathway in comparison to the conventional approach based on the pore geometry. For small gas molecules, the potential energy barrier for diffusion through a nanostructured material is not due to the excluded volume effects as often assumed in conventional methods but primarily due to its strong binding with the porous material. As a result, the minimum energy diffusion pathway must be identified from the potential energy gradient. Toward that end, the string method is proved to be both more efficient and accurate than the geometric analysis. We have shown, with the help of the string method, that TST is able to predict the selfdiffusion coefficients of methane and hydrogen in a number of MOFs and zeolites in excellent agreement with results from direct MD simulation. By contrast, good agreement was

13 ACS Paragon Plus Environment

Langmuir

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 24

accomplished in previous work by applying a dynamic correction factor9, 11. While dynamic correction is inevitable for understanding diffusion of larger molecules at finite loadings, it seems that its effect is relatively insignificant for small gas molecules at infinite dilution. Although its extension to non-spherical gases in nanostructure materials is yet to be established, classical density functional theory (cDFT) has been proved efficient for calculating the adsorption isotherms and thermodynamic properties of methane, hydrogen and carbon dioxide in a large library of nanostructured materials31. In combination with the transition state theory (TST) and excess entropy scaling, we demonstrate that the thermodynamic route is efficient for predicting the self-diffusion coefficients at both infinite dilution and finite loadings. The new procedure is able to improve the computational efficiency by few orders of magnitude in comparison with MD simulation. Whereas the excess entropy scaling method was reported in our earlier publication, accurate evaluation of gas diffusivity at infinite dilution is computationally challenging in its own right and, as shown in this work, the string method obviates some of the main difficulties in application of the transition state theory. In this work, we consider self-diffusivity in nanoporous materials at conditions commonly used for gas storage and separation. However, the theoretical procedure is expected broadly applicable to arbitrary conditions including systems at low temperature where gas molecules exhibit capillary condensation. Because the string method and the transition state theory are concerned only with a single gas molecule in a specific nanoporous material, computation of gas diffusivity at infinite dilution is independent of the phase behavior of the adsorbate. Furthermore, the excess entropy scaling method is valid for both gas and liquid systems33, 34. Given the feasibility of the cDFT method to predict thermodynamic properties and

14 ACS Paragon Plus Environment

Page 15 of 24

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

phase behavior at both supercritical and subcritical conditions, we anticipate that its combination with the excess entropy scaling will yield reliable diffusivity at finite loading. Many nanoporous materials can be synthesized via modular construction from virtually infinite types of linkers and nodes. Due to the vast number of possible structures and their underlying potential for versatile applications in gas storage and separation, it is crucial to develop more efficient computational methods to predict and compare the thermodynamic and transport properties of gas molecules in different porous materials. Toward that end, it is our hope that the theoretical procedure introduced in this work will help establish structure-property relationships useful for computational material design, synthesis as well as process optimization. Acknowledgments The authors are grateful to the U.S. National Science Foundation (NSF-CBET-1404046) for the financial support. The simulation work was performed at the National Energy Research Scientific Computing Center (NERSC), supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC02-05CH11231. References 1. Karger, J.; Ruthven, D. M.; Theodorou, D. N., Diffusion in Nanoporous Materials. Wiley-VCH: Weinheim, Germany, 2012. 2. Krishna, R., Diffusion in porous crystalline materials. Chem. Soc. Rev. 2012, 41, (8), 3099-3118. 3. Smit, B.; Maesen, T. L. M., Molecular Simulations of Zeolites: Adsorption, Diffusion, and Shape Selectivity. Chem. Rev. 2008, 108, (10), 4125-4184. 4. Karger, J.; Ruthven, D. M., Diffusion in nanoporous materials: fundamental principles, insights and challenges. New J. Chem. 2016, 40, (5), 4027-4048. 5. Dubbeldam, D.; Snurr, R. Q., Recent developments in the molecular modeling of diffusion in nanoporous materials. Molecular Simulation 2007, 33, (4-5), 305-325. 6. Skoulidas, A. I.; Sholl, D. S., Molecular Dynamics Simulations of Self-Diffusivities, Corrected Diffusivities, and Transport Diffusivities of Light Gases in Four Silica Zeolites To Assess Influences of Pore Shape and Connectivity. J. Phys. Chem. A 2003, 107, (47), 1013210141. 7. Auerbach, S. M., Theory and simulation of jump dynamics, diffusion and phase equilibrium in nanopores. Int. Rev. Phys. Chem. 2000, 19, (2), 155-198. 15 ACS Paragon Plus Environment

Langmuir

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 24

8. Beerdsen, E.; Dubbeldam, D.; Smit, B., Molecular Understanding of Diffusion in Confinement. Phys. Rev. Lett. 2005, 95, (16), 164505. 9. Dubbeldam, D.; Beerdsen, E.; Vlugt, T. J. H.; Smit, B., Molecular simulation of loadingdependent diffusion in nanoporous materials using extended dynamically corrected transition state theory. J. Chem. Phys. 2005, 122, (22), 224712. 10. Haldoupis, E.; Nair, S.; Sholl, D. S., Efficient Calculation of Diffusion Limitations in Metal Organic Framework Materials: A Tool for Identifying Materials for Kinetic Separations. J. Am. Chem. Soc. 2010, 132, (21), 7528-7539. 11. June, R. L.; Bell, A. T.; Theodorou, D. N., Transition-state studies of xenon and sulfur hexafluoride diffusion in silicalite. J. Phys. Chem. 1991, 95, (22), 8866-8878. 12. Tunca, C.; Ford, D. M., A transition-state theory approach to adsorbate dynamics at arbitrary loadings. J. Chem. Phys. 1999, 111, (6), 2751-2760. 13. Beerdsen, E.; Smit, B.; Dubbeldam, D., Molecular simulation of loading dependent slow diffusion in confined systems. Phys. Rev. Lett. 2004, 93, (24), Artn 248301. 14. Camp, J. S.; Sholl, D. S., Transition State Theory Methods To Measure Diffusion in Flexible Nanoporous Materials: Application to a Porous Organic Cage Crystal. J. Phys. Chem. C 2016, 120, (2), 1110-1120. 15. Beerdsen, E.; Dubbeldam, D.; Smit, B., Loading Dependence of the Diffusion Coefficient of Methane in Nanoporous Materials. The Journal of Physical Chemistry B 2006, 110, (45), 22754-22772. 16. Ruiz-Montero, M. J.; Frenkel, D.; Brey, J. J., Efficient schemes to compute diffusive barrier crossing rates. Mol. Phys. 1997, 90, (6), 925-941. 17. E, W.; Ren, W. Q.; Vanden-Eijnden, E., String method for the study of rare events. Phys. Rev. B 2002, 66, (5), Artn 052301. 18. Kolokathis, P. D.; Káli, G.; Jobic, H.; Theodorou, D. N., Diffusion of Aromatics in Silicalite-1: Experimental and Theoretical Evidence of Entropic Barriers. J. Phys. Chem. C 2016, 120, (38), 21410-21426. 19. Rosenfeld, Y., Relation between Transport-Coefficients and Internal Entropy of Simple Systems. Phys. Rev. A 1977, 15, (6), 2545-2549. 20. Liu, Y.; Fu, J.; Wu, J., Excess-Entropy Scaling for Gas Diffusivity in Nanoporous Materials. Langmuir 2013, 29, (42), 12997-13002. 21. Eddaoudi, M.; Kim, J.; Rosi, N.; Vodak, D.; Wachter, J.; O'Keeffe, M.; Yaghi, O. M., Systematic Design of Pore Size and Functionality in Isoreticular MOFs and Their Application in Methane Storage. Science 2002, 295, (5554), 469-472. 22. Chui, S. S.-Y.; Lo, S. M.-F.; Charmant, J. P. H.; Orpen, A. G.; Williams, I. D., A Chemically Functionalizable Nanoporous Material [Cu3(TMA)2(H2O)3]n. Science 1999, 283, (5405), 1148-1150. 23. Rappe, A. K.; Casewit, C. J.; Colwell, K. S.; Goddard, W. A.; Skiff, W. M., UFF, a full periodic table force field for molecular mechanics and molecular dynamics simulations. J. Am. Chem. Soc. 1992, 114, (25), 10024-10035. 24. Xu, X.; Ting, C. L.; Kusaka, I.; Wang, Z.-G., Nucleation in Polymers and Soft Matter. Annu. Rev. Phys. Chem. 2014, 65, (1), 449-475. 25. E, W.; Ren, W.; Vanden-Eijnden, E., Simplified and improved string method for computing the minimum energy paths in barrier-crossing events. The Journal of Chemical Physics 2007, 126, (16), 164103.

16 ACS Paragon Plus Environment

Page 17 of 24

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

26. Fu, J.; Liu, Y.; Tian, Y.; Wu, J. Z., Density Functional Methods for Fast Screening of Metal Organic Frameworks for Hydrogen Storage. J. Phys. Chem. C 2015, 119, (10), 5374-5385. 27. Haldoupis, E.; Nair, S.; Sholl, D. S., Pore size analysis of >250 000 hypothetical zeolites. PCCP 2011, 13, (11), 5053-5060. 28. Frenkel, D.; Smit, B., Understanding Molecular Simulation: from Algorithms to Applications. 2nd ed.; Academic Press: San Diego, 2002; p 638. 29. Kolokathis, P. D.; Pantatosaki, E.; Gatsiou, C.-A.; Jobic, H.; Papadopoulos, G. K.; Theodorou, D. N., Dimensionality reduction of free energy profiles of benzene in silicalite-1: calculation of diffusion coefficients using transition state theory. Molecular Simulation 2014, 40, (1-3), 80-100. 30. Ford, D. C.; Dubbeldam, D.; Snurr, R. Q., The effect of framework flexibility on diffusion of small molecules in the metal-organic framework IRMOF-1. Diffusionfundamentals.org 2009, 11, (78), 1-8. 31. Fu, J.; Tian, Y.; Wu, J., Classical density functional theory for methane adsorption in metal-organic framework materials. AlChE J. 2015, 61, (9), 3012-3021. 32. Liu, B.; Yang, Q.; Xue, C.; Zhong, C.; Smit, B., Molecular simulation of hydrogen diffusion in interpenetrated metal-organic frameworks. PCCP 2008, 10, (22), 3244-3249. 33. Hopp, M.; Gross, J., Thermal Conductivity of Real Substances from Excess Entropy Scaling Using PCP-SAFT. Ind. Eng. Chem. Res. 2017, 56, (15), 4527-4538. 34. Sahu, P.; Ali, S. M.; Shenoy, K. T., Test of Universal Scaling Law for Molecular Diffusion of Liquids in Bulk and Nanotube Confinement. J. Phys. Chem. C 2017, 121, (21), 11968-11974.

17 ACS Paragon Plus Environment

Langmuir

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure Captions Figure 1. a) Reaction coordinate diagrams for CH4 diffusion along the major pore of MOF-5 calculated from the geometric analysis (dashed line) and from the string method (solid line). The numbered circles represent the minimum and maximum energy states along the pore center (open circles) and the diffusion pathway (filled circles). b) The locations of CH4 molecule within the 3-dimensional structure of MOF-5 corresponding to the extremum energy states. The framework atoms are shown as small spheres: Zn, purple; O, red; C, gray; H, white.

Figure 2. Self-diffusivity of H2 and CH4 at 298 K and infinite dilution in several MOFs. The black squares are TST predictions from this work, the MD simulation data and results from a previous implementation of TST are from the literature10.

Figure 3. Self-diffusivity of CH4 in six zeolites at 298 K and infinite dilution from the TSTstring method (black squares), MD simulation (dashed line), and previous TST predictions (red spheres)27. The 3-letter codes stand for different types of zeolites.

Figure 4. Comparison of the excess entropy scaling and MD simulation for the self-diffusivity of CH4 at 298 K in MOF-5 (a) and in CuBTC (b). The MD-DFT predictions correspond to excess entropy scaling based on the MD simulation for diffusivity at infinite dilution. The MD results are from the literature6.

Figure 5. Adsorption isotherms for CH4 in MOF-5 (a) and in CuBTC (b) at 298K. The lines are from cDFT predictions, and squares are from GCMC simulation31.

18 ACS Paragon Plus Environment

Page 18 of 24

Page 19 of 24

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

Figure 1

19 ACS Paragon Plus Environment

Langmuir

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 2

20 ACS Paragon Plus Environment

Page 20 of 24

Page 21 of 24

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

Figure 3

21 ACS Paragon Plus Environment

Langmuir

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 4

22 ACS Paragon Plus Environment

Page 22 of 24

Page 23 of 24

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

Figure 5

23 ACS Paragon Plus Environment

Langmuir

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19

ACS Paragon Plus Environment

Page 24 of 24