Design Principles for the Ultimate Gas Deliverable Capacity Material: Nonporous to Porous Deformations Without Volume Change

Understanding the fundamental limits of gas deliverable capacity in porous materials is of critical importance as it informs whether technical targets (e.g., for on-board vehicular storage) are feasible. High-throughput screening studies of rigid materials, for example, have shown they are not able to achieve the original ARPA-E methane storage targets, yet an interesting question remains: what is the upper limit of deliverable capacity in flexible materials? In this work we develop a statistical adsorption model that specifically probes the limit of deliverable capacity in intrinsically flexible materials. The resulting adsorption thermodynamics indicate that a perfectly engineered, intrinsically flexible nanoporous material could achieve higher methane deliverable capacity than the best benchmark systems known to date with little to no total volume change. Density functional theory and grand canonical Monte Carlo simulations identify a known metal-organic framework (MOF) that validates key features of the model. Therefore, this work (1) motivates a continued, extensive effort to rationally design a porous material analogous to the adsorption model and (2) calls for continued discovery of additional high deliverable capacity materials that remain hidden from rigid structure screening studies due to nominal non-porosity. Abstract Understanding the fundamental limits of gas deliverable capacity in porous materials is of critical importance as it informs whether technical targets (e.g., for on-board vehicular storage) are feasible. High-throughput screening studies of rigid materials, for example, have shown they are not able to achieve the original ARPA-E methane storage targets, yet an interesting question remains: what is the upper limit of deliverable capacity in flexible materials? In this work we develop a statistical adsorption model that specifically probes the limit of deliverable capacity in intrinsically flexible materials. The resulting adsorption thermodynamics indicate that a perfectly en-gineered, intrinsically flexible nanoporous material could achieve higher methane deliverable capacity than the best benchmark systems known to date with little to no total volume change. Density functional theory and grand canonical Monte Carlo simulations identify a known metal-organic framework (MOF) that validates key features of the model. Therefore, this work (1) motivates a continued, extensive effort to rationally design a porous material analogous to the adsorption model and (2) calls for continued discovery of additional high deliverable capacity materials that remain hidden from rigid structure screening studies due to nominal non-porosity.

Understanding the fundamental limits of gas deliverable capacity in porous materials is of critical importance as it informs whether technical targets (e.g., for on-board vehicular storage) are feasible. High-throughput screening studies of rigid materials, for example, have shown they are not able to achieve the original ARPA-E methane storage targets, yet an interesting question remains: what is the upper limit of deliverable capacity in flexible materials? In this work we develop a statistical adsorption model that specifically probes the limit of deliverable capacity in intrinsically flexible materials. The resulting adsorption thermodynamics indicate that a perfectly engineered, intrinsically flexible nanoporous material could achieve higher methane deliverable capacity than the best benchmark systems known to date with little to no total volume change. Density functional theory and grand canonical Monte Carlo simulations identify a known metal-organic framework (MOF) that validates key features of the model. Therefore, this work (1) motivates a continued, extensive effort to rationally design a porous material analogous to the adsorption model and (2) calls for continued discovery of additional high deliverable capacity materials that remain hidden from rigid structure screening studies due to nominal non-porosity.

Introduction
The economic viability of gas storage technologies for cleaner vehicle transportation hinges on materials with sufficiently large deliverable capacity, defined as the difference in gas loading between the charging and discharging pressures. [1][2][3][4] There have therefore been concerted efforts to understand the physical limit of deliverable capacity in porous materials, especially within the applications of hydrogen and methane storage. [5][6][7][8] We will provide evidence that, through a carefully designed analytical model and corroborating calculations on known materials, the deliverable capacities of today's benchmark materials still have not reached a fundamental limit if intrinsic flexibility can be deftly exploited. While this work uses methane storage as a case study, the lessons learned are general and could be applied to other gas stor-age applications. 9 Typical computational exploration of gas deliverable capacity in porous materials consists of screening large databases of known and hypothetical materials using the rigid structure approximation. [5][6][7][8] For example, these studies have convincingly revealed that the original Advanced Research Project Agency-Energy (ARPA-E) methane deliverable capacity (DC) targets of 315 v(STP)/v delivered at ambient temperature between 5.8 and 65 bar are not achievable using the materials which have been screened under the rigid framework approximation. 7 One fundamental problem is that strong methane interactions are needed to achieve high capacities at the charging pressure, but this results in a capacity that is also too high at the discharging pressure, thereby reducing the deliverable capacity to ∼ 70% or less of the high pressure adsorption. This realization has also been confirmed by a host of modelling studies that elucidate the deliverable capacity upper bound in rigid materials. 10,11 Temperature swing operation can address this decrease in deliverable capacity, 12 but at the cost of additional weight for the necessary heating/cooling systems. Extrinsically flexible materials on the other hand, 13,14 which deliver gas via an S-shaped isotherm (e.g., IUPAC Type V isotherm) due to a large volume expansion from a nonporous to a porous state, have been investigated thoroughly both computationally and experimentally in an attempt to overcome this physical limit. [15][16][17][18] However, it is unclear how well these large volume change materials can be practically utilized from a systems design standpoint, as they must remain functional after many cycles despite shear displacive phase transformations. 19,20 Current technology development efforts appear to forsake the enhanced deliverable capacity of extrinsically flexible materials, and instead focus on constant volume, rigid adsorbents that can be densified with robust mechanical stability. 21 In the search for materials that improve upon state-of-the-art methane DC targets, here we diverge from the approaches above.
First, we design an adsorption model of an intrinsically flexible material 22-25 (exhibits a pore volume change but no significant unit cell volume change) that offers the same favorable adsorption properties of an extrinsically flexible material. We describe a flexible slit pore model that provides a physical basis for designing high deliverable capacity materials, whereby an "ideal" material could translate nearly all of the adsorption capacity to deliverable capacity (see Figure 1 for an example). We specifically design the minimum energy configuration of the pore to be too small to accommodate adsorbates, and the fully porous state is accessible at higher energy, which can be stabilized by the presence of adsorbates. 26 The deliverable capacity can therefore be enhanced due to the reduced uptake at low pressures (while maintaining an equally high uptake at high pressures) to the point that it competes with or even exceeds today's benchmark materials. 18 Figure 1: Example rigid material isotherms where adsorption is either too energetically favorable (blue), optimized (orange), or too unfavorable (green) for maximizing DC. Flexible systems (black) have the potential to maximize DC by suppressing adsorption at low pressures without sacrificing high capacity at high pressures.
Second, rather than computationally screening hundreds of thousands of materials within the rigid structure approximation, we perform a set of detailed density functional theory (DFT) and grand canonical Monte Carlo (GCMC) calculations on a family of intrinsically flexible materials, known as M(NDC) (M={Ca,Sr}, NDC=1,4-naphthalenedicarboxylate). 28 DFT calculations validate some of the key geometric/energetic features of the flexible slit pore model that impart its high deliverable capacity.
We initially selected this material because DFT relaxation of the desolvated structure resulted in a large pore volume change but very little unit cell volume change. Despite their "nominal non-porosity" in the minimum energy state (i.e. the pore size of the minimum energy configuration of the framework is too small to accommodate methane), these materials are subject to significant changes in the pore size/shape with minimal change to the unit cell volume. The transition from nonporous to porous becomes energetically favorable due to adsorption. 25 In similar systems, cooperative adsorption due to ring rotations in constant volume materials has been modelled and observed (although the minimum energy state was still porous). 29,30 It should be noted that such materials could only ever be computationally identified with a fully flexible treatment coupling framework dynamics and adsorption; 16,31-34 hence, they will have remained hidden from almost any high-throughput adsorbent screening study to date. 35,36 Finally, we show that the isoreticular tunability 37,38 of metal-organic frameworks (MOFs) can be used to design an improved version of the M(NDC) adsorbent. Ultimately, this work outlines a path forward to discover deliverable capacity materials through concrete design strategies motivated by the models and simulations developed herein.

Flexible Slit Pore Model
We derive a numerically-solvable adsorption model that accounts for intrinsic flexibility via a rotating and vibrating slit pore. This allows the adsorption properties to be quickly and systematically solved over the model parameter space and provides insights into deliverable capacity limits in intrinsically flexible materials. Figure 2 visualizes our slit pore adsorption model where the centerlines of the pore walls are separated by a distance L. The tilt of the pore, θ, is the angle between the centerline and the normal to the i crystallographic vector. x, the adsorbate's degree of freedom, is the magnitude of its displacement vector from the centerline of the left hand wall in Figure 2, parallel with i. The slit pore is large enough to ignore edge effects such that the particle's adsorption energy will be independent of its coordinates in the j and k (into the page) directions. x ⊥ is the magnitude of the displacement vector normal to the centerline of the wall. Therefore the minimum distance to the first wall is given by x ⊥,1 (x, θ) = x cos θ and the minimum distance to the opposite wall by x ⊥,2 (x, θ, L) = (L − x) cos θ. The largest included sphere that can fit without overlapping the pore wall is D i (θ, L) = L cos θ − σ F , where σ F is the Lennard-Jones diameter of the pore wall.

Slit pore geometry
This model will enable prediction of adsorption thermodynamics with numerical integration rather than Monte Carlo simulations, and despite its simplicity, an exact analog of this rotating slit pore can be found within a known MOF (visualized in Figure 2 and discussed later). At θ = 0, L can be chosen such that D i is slightly larger than the adsorbate and maximizes the binding strength. Then, if L remains constant, increasing θ decreases D i until the slit pore becomes nonporous. If the nonporous state is energetically preferred but the porous state can be stabilized by the presence of adsorbates, a significant pore volume change can be achieved without a significant cell volume change. This has the potential to provide enhanced DC, which we quantitatively explore with the model.

Framework-guest potential
Lennard-Jones interactions describe frameworkguest interactions, and the methane guest parameters ϵ G = 148 K and σ G = 3.73 Å are taken from TraPPE. 39 The framework consists of carbon atoms with ϵ F = 52.8 K and σ F = 3.43 Å taken from the Universal Force Field. 40 Each pore wall is treated as a plane of uniformly smeared Lennard-Jones spheres with planar density, d s , equal to graphene. The potential energy between methane and one of the smeared carbon walls, E s , can then be derived  Figure 2: The relevant geometry parameters for computing the energetics between an adsorbate (gray circle) and the slit pore walls (green rectangles). Red dashed lines represent the center line of each wall, red dots represent the rotation center for each wall, and the black dots represents the adsorbate center. We also show example of a known MOF, Sr(NDC), we discovered to have rotating slit pores.
as a function of x ⊥ , 41 .

Framework potential
We consider three different forms of the framework potential energy, E F , corresponding to the different types of flexibility given by eqs. (3) to (5).
1. In the static (S) slit pore of eq. (3), θ and L are fixed at θ eq and L eq and do not fluctuate.
2. In the rotating (R) slit pore of eq. (4), θ fluctuates in a skewed double well potential, the shape of which is controlled by V max , a, b, and c. L is fixed at L eq and does not a fluctuate.
3. In the rotating and vibrating (R+V) slit pore of eq. (5), L fluctuates according to a harmonic spring potential with spring constant k.
Choosing V max = 20 kJ/mol, a = π/6, b = π/12, c = 50 kJ/mol creates an E F profile similar to the MIL-53 breathing free energy profile, 42 and the effect of parameter variation on the adsorption performance is examined further in S1.1. The global minimum of E R F with these parameters occurs at θ ′ = 31.1 • and corresponds to an asymmetric double well, 43 although different parameter sets can yield asymmetric single well potentials and still produce similar favorable adsorption properties (see S1.1).

Adsorption volume discretization
Next we discretize the adsorption space into lattice sites and disallow multi-site occupancy to model repulsive particle overlaps. Note that site discretization is similarly used to derive the Langmuir isotherm from statistical mechanics. 10 Specifically, the adsorption space is discretized into M channels in the j,k plane as shown in Figure 3. Each channel site has = !"#$ = % Figure 3: Discretization of the adsorption domain from Figure 2 into M = 4 lattice channels in the (j,k) plane with N = 2 adsorbates, each with total and free areas, A site and A o , respectively. a cross-sectional area of A site = (2 1/6 σ G ) 2 = (4.19 Å) 2 with a volume of V site = A site L. The calculation of the free cross-sectional area, A o = 1.1 Å 2 , is outlined in S1.3 and its accuracy is ensured later (see Figure 6). While preventing particle overlaps, this discretization does not account for favorable guest-guest interactions. However, the adsorbates in the slit pore have a limited number of nearest neighbors and each methane-methane interaction could contribute at most ϵ = −1.2 kJ/mol (minimum in the methane-methane interaction potential), making it a small contribution relative to the interaction between methane and the framework in this particular system. Guest-guest interactions could of course be accounted for explicitly, but then we would rely on MC simulations to solve the model due to the complexity of the resulting partition function. Therefore, the single site partition function for an adsorbed methane molecule in one of the site channels becomes where Λ is the thermal de Broglie wavelength.

Canonical partition function
The canonical partition functions corresponding to each type of flexible slit pore are expressed as a product of single site partition functions with integration over the framework degrees of freedom, θ and L, in eqs. (7) to (9).
Here β = (k B T ) −1 where k B is Boltzmann's constant and T is temperature, and ( M N ) accounts for the indistinguishability of particles and enforces single occupancy sites. These one, two, or three dimensional integrals can be evaluated using standard numerical integration libraries contained in scipy or Mathematica (see S1.4 for code availability).

Isotherms, adsorption enthalpy, and deliverable capacity
Solving the canonical partition function permits calculation of the grand partition function. Then the adsorption isotherm, or the expectation number of particles, ⟨N ⟩, as a function of µ, the chemical potential, can be computed as A material's DC is the difference in adsorbate loading between the adsorption pressure, P ads , and the desorption pressure, P des < P ads , We relate µ to P through the ideal gas law, 44 a good approximation at ambient temperatures for light gases such as methane. For the DOE storage targets, P des = 5.8 bar and P ads = 65 bar (see Figure 1). Finally, we can compute the differential enthalpy of adsorption via the Clausius-Clapeyron formalism from isotherms obtained at multiple temperatures, where R is the gas constant and P 0 = 1 bar is the ideal gas reference state pressure. ∆h gives the amount of heat that must be added/removed from the system to maintain constant temperature when transferring a molecule from the gas phase to the adsorbent, or vice versa, at a given loading. 45,46

DFT flexibility calculations
Rigorous atomistic simulations of adsorption in flexible nanoporous materials are complicated endeavors, not only because they often require enhanced sampling techniques, 32,34 but also because the predictions can be very sensitive to small fluctuations in structural dynamics. 35 In other words, quantitative isotherm predictions require an extremely accurate description of the potential energy surface and long simulations. Each challenges either classical force fields (do not universally transfer with high accuracy across the diverse MOF chemistry/structure space) or ab initio molecular dynamics (are very expensive). We therefore turn to a simple DFT screening procedure to gauge if a material with a rotating slit pore has potential for an intrinsic nonporous to porous transitions. We calculate two contributions to the potential energy in a nonporous to porous transition, the framework deformation energy and the adsorption energy, by considering the empty framework, F, or the framework with guest adsorbates, F · G. The system of interest X is labelled by X| Z Y , where Y is the system from which the geometry of X was taken. Z ∈ {DFT, Exp.} denotes whether the geometry of Y was determined by a DFT relaxation or taken from experimental measurements. So, F| DFT F·G denotes that we perform a DFT optimization of the framework with adsorbates present but then keep only the framework. F| Exp.
F·G denotes we took the experimental CIF file and removed the adsorbates/solvent that were present, and F| DFT F indicates we have DFT optimized the empty framework. The DFT computed energy of these configurations is denoted E(X| Z Y ), and so the total energy required to deform the framework to accommodate a given set of adsorbates, ∆E def , the total binding energy of those adsorbates, ∆E b , and the net energy change ∆E net are expressed as These can be normalized by the N number of adsorbates in the F · G configuration:

Results
First, we predict the adsorption properties of the flexible slit pore model. Mapping the model's DC to a volumetric basis then provides useful insights into the maximum achievable DC in intrinsically flexible materials and whether such a material could exceed the performance of the current benchmark systems. Next, we utilize DFT to probe the adsorption energy landscape of a known material that exhibits some of the advantages featured by the adsorption model. Finally, we perform an isoreticular design exercise to improve its potential for methane storage.

Adsorption in flexible slit pores
The main limitation of rigid materials 6,7,10 is briefly summarized by the static pore isotherms in Figure 4, where the deliverable capacity is calculated as a function of pore wall separation distance (at fixed θ eq = 0). Two maximums in the deliverable capacity occur when the methane binding energy is neither too strong nor too weak. In the intermediary local minimum, the adsorption is maximized at 65 bar due to highly favorable methane interactions, but the DC suffers because too much methane is also adsorbed at 5.8 bar. This represents the fundamental limitation of gas storage in rigid materials, and, even for the record methane storage in rigid materials, the deliverable capacity is ∼ 30% less than the adsorption capacity at 65 bar. 21 In contrast, the R and R+V slit pores avoid this trade-off and allow nearly 100% of the material's adsorption capacity to be utilized as deliverable capacity. When L eq is too small, even a pore tilt of θ = 0 yields a nonporous material, and so the uptake and DC are negligible. When L eq is too large, we recover Langmuirian adsorption because the framework is porous for any value of θ. The "sweet-spot" for achieving an optimal DC is a mere ∼0.5 Å window in L eq . At these separation distances, the framework adopts a nonporous state at θ ′ in the absence of any adsorbates, while the porous state (θ ∼ 0) is stabilized at sufficiently high pressures by the presence of adsorbates. This leads to an inflection in the isotherm, whereby adsorption is suppressed at low pressures and deliverable capacity is increased. A critical note is that this S-shaped isotherm behavior can only be achieved in this system due to "cooperative adsorption", i.e. if M > 1. In S1.2 we elaborate more on this topic and demonstrate how DC varies with M , i.e. the number of adsorption sites that exist in the open slit pore. In other statistical mechanics models of flexible adsorbents, this cooperativity manifests through different mechanisms but still causes inflections in the isotherm. 26,29 Framework vibrations at finite temperatures can have a non-negligible impact on the pore size distribution and therefore a porous material's separation capabilities, 35 so it is important to note that the optimal deliverable capacity of the R+V pore is unaffected, albeit at occurring at slightly lower L eq .

Intrinsic heat management
The flexible slit pore model demonstrates intrinsic heat management in the same way as extrinsically flexible materials. 18 The positive energy penalty required to deform the adsorbent offsets some of the energy released during an adsorption event, leading to a reduction in the differential enthalpy of adsorption. Figure 5 shows ∆h computed for the R and S slit pores. In this plot, θ eq = 0 for S. For both S and R, the L eq that maximizes deliverable capacity is used (6.7 and 7.2 Å respectively from Figure 4). Since the static slit pore corresponds to Langmuirian adsorption, ∆h is constant as a function of N . While the rotating slit pore still shares some Langmuirian assumptions, the underlying E F profile results in a reduced ∆h. This reduced ∆h results in less heat that needs to be added/removed during desorption/adsorption cycling, thereby re- ducing space, cost, and engineering complexity associated with on-board thermal management systems.

Volumetric deliverable capacity
The technical targets for gas deliverable capacity are often specified on a volumetric basis, i.e., the original ARPA-E targets required a deliverable capacity of 315 v(STP)/v. Therefore we convert the model's adsorbate density from a molecules per site basis to a volumetric basis via where P atm and T are standard pressure and temperature and N av is Avogadro's number. We obtain saturation capacity on a volumetric basis, σ vol , if we take ⟨N ⟩/M → 1 in eq. (14), which corresponds to taking P ads → ∞. Replotting the volumetric adsorbate densities and volumetric deliverable capacity of the static and rotating slit pores in Figure 6 reveals several important insights. First, the charging pressure capacity of our idealized static slit pore model, ρ vol at P ads , is consistent with GCMC simulations of methane adsorption in a bilayer graphene structure as a function of separation distance, L eq , between each sheet (see S2.2). This indicates that the model's site volume discretization closely reproduces the adsorption of an analogous atomistic simulation. Second, the static pore exhibits an optimal DC vol ∼ 170 v(STP)/v, and the highest reported deliverable capacities for rigd MOFs are ∼200 v(STP)/v, again validating the adsorption volume discretization. In the case of the top-performing HKUST-1, only 75 % of the 65 bar adsorption capacity, ρ vol ∼270 v(STP)/v, is converted to deliverable capacity, DC vol ∼200 v(STP)/v; and these numbers are only further reduced when MOF densification is required. Most importantly, since the rotating slit pore model can translate ∼95 % of its adsorption capacity to DC (Figure 4), a porous material that exactly mimics the geometry and potential energy surface of the model could significantly outperform the best known rigid materials. The caveat is of course that the model has a highly idealized geometry that will be difficult to engineer exactly into a synthesizable crystal structure, and DC vol is still constrained by the model's upper bound on adsorption capacity of ∼290 v(STP)/v.

Nonporous to porous without volume change in M(NDC)
The  Figure 7). If we instead replace each DMF with a methane molecule, DFT relaxation maintains the same porous state where all slit pores remain aligned in parallel. Deleting the DMF from the experimental structure (F| Exp. F·8DMF ) or deleting methane from the DFT optimized structure (F| DFT F·8CH 4 ) yields pores that can accommodate methane, while F| DFT F is nonporous. Table 1 shows the Zeo++ computed geometric properties of these configurations where D i is the largest included sphere, D f is the largest free sphere, V is the crystal cell volume, V 0 is the crystal cell volume of the F| DFT F configuration, and V occ /V is "probe occupiable" (probe diameter = 3.4 Å) pore vol-ume 58,59 fraction that includes both the accessible and non-accessible volumes. Notably, V /V 0 expansion is practically negligible in comparison to other phase change methane adsorbents, whereas V occ /V increases drastically. The de-

Config
.  To gauge whether these materials have favorable methane adsorption properties requires more than the geometric analysis above. Namely, we compute energetics of the adsorbate-driven porous to nonporous transition described by (13) and summarize them in Table 2. Several F·G configurations are considered, with G ∈ {1CH 4 , 8CH 4 , 16CH 4 }. In 1CH 4 , the structure is optimized with just one methane per simulation box to probe the binding energy and framework deformation energy associated with the first adsorption event. 8CH 4 corresponds to one methane at each open metal site, and 16CH 4 corresponds to an overall density of two methane molecules per open metal site, i.e. saturation loading. Despite the energy penalty required to deform the framework from the nonporous to porous state, the net adsorption process is favorable for all loadings considered. Importantly, a "collective" behavior is evident Table 2. ∆Ē def is large for G = 1CH 4 since it induces a local deformation, i.e., opening of a single slit pore. This energy penalty drops significantly when all slit pores are aligned in the open con-  16CH 4 ), albeit at slightly lower binding energy per methane. The final result is that the net energy decrease of adsorbing methane is more favorable near saturation than in the empty framework. The ratio of cell volumes between the porous and nonporous states, V /V o , is very close to 1. However, the difference in GCMC predicted isotherms between the F| DFT F and F| DFT F·8CH 4 configurations in Sr(NDC) is drastic, as shown in Figure 8, due to the large increase in V occ /V . Force field information and GCMC settings are provided in S2.2. 39,40,60 How, and at what pressure, any S-shaped adsorption profile occurs will clearly depend on the underlying framework and framework-guest potential energy surface, although we are unaware of any specialized MOF force fields 61-64 designed specifically for this Sr coordination environment that also capture the non-bonded interactions between methane and the open metal site. [65][66][67][68] Nonetheless, the true isotherm will be bounded between the F| DFT F and F| DFT F·8CH 4 isotherms. Therefore, the deliverable capacity cannot exceed ∼150 v(STP)/v.

Experimental Activation
Sr(NDC) was readily synthesized in large quantities (see S3) and good agreement was obtained between the powder x-ray diffraction (PXRD) of the as-synthesized sample and the previously reported structure (see Figure 9 and S3 for experimental details). Extreme activation condi-  tions were originally found to be necessary to overcome the strong binding between the DMF solvent and the open metal site and allow it to release through the narrow pore constrictions in the framework, 28 and we found that heating at 300 • C for 16 hours led to a structural change in the PXRD pattern of the activated sample. The primary reflection is preserved, suggesting that the unit cell remains unchanged. Nonetheless, while Sr(NDC) may be an academically interesting manifestation of the slit pore model, more work remains to discover a slit pore material whose D f is not prohibitively small for methane adsorption.

Isoreticular expansion of M(NDC)
The deficiencies with Sr(NDC), and where it differs from the adsorption model, are twofold. First, the individual slit pores are staggered slightly, which leads to a D i ̸ = D f . This problematically low D f means that the barrier to methane diffusion will be large and points to kinetics-limited adsorption (note this manifests in the extremely high desorption temperatures, ∼ 300 • C, required to remove DMF and activate the material). 28 The second issue is that the probe accessible pore volume of the open structure is not large enough, evidenced by the ∼150 v(STP)/v saturation ca-  pacity of the F| DFT F·G isotherm in Figure 8. More void volume could be imparted by an isoreticular expansion with a larger planar linker, as long as the energy minimized state in the absence of adsorbates remains nonporous. We performed isoreticular expansion of Sr(NDC) using hypothetical linkers and a 1D rod MOF assembly algorithm, 69 then repeated our porosity and adsorption energetics analysis (detailed in S2.3). Even though the isoreticular analogs indeed show the same slit pore rotation, the expanded linkers introduce varying degrees of porosity into the F| DFT F state, which limits DC improvements. Nonetheless, it demonstrates that careful rational design of MOF structures can be useful, and perhaps will even be necessary, to design the ultimate deliverable capacity materials via this flexible slit pore approach.

Discussion and Conclusions
The results of our intrinsically flexible slit pore model highlight some interesting challenges and opportunities in the continual search for adsorbents with optimal gas deliverable capacity. The deliverable capacities of the rotating slit pore model are, at minimum, a promising indicator that the best gas storage materials may be (as of yet) undiscovered and surpass benchmark rigid and extrinsically flexible adsorbents. This arises from an adsorption induced, nonporous to porous transition that can occur little to no unit cell volume change. The first challenge in accomplishing this is alluded to in Figure 4. The L eq window over which the optimal rotating slit pore exists is ∼0.5 Å. In other words, a computational approach to quantitatively screen potential adsorbents must have an extremely accurate potential energy surface. While significant progress has been made in MOF force field development, no single universal force field exists that is highly accurate for all the complex chemistry and topologies exhibited by MOFs. Second, new materials would have to be tested or simulated that are nominally nonporous. Fully flexible simulations are required to evaluate their potential, 33,34 and therefore would have remained undiscovered from rigid screening studies, for which there are very few exceptions. 35,36 Third, one must find an adsorbent that mimics the idealized flexible slit pore model as closely as possible in order to maximize deliverable capacity.
These challenges notwithstanding, we utilized the insights from the adsorption model to identify a known MOF that displays most, but not all, of the advantages of the slit pore model. Namely, DFT calculations reveal it exhibits a nonporous minimum energy state, and undergoes an adsorbate-driven, energetically favorable, and constant volume transformation to a porous state due to a slit pore rotation. We additionally noted that both D f and V occ /V are too small and constitute sufficient limitations that Sr(NDC) cannot achieve a record deliverable capacity. We sought to address these shortcomings via isoreticular expansion of the framework with a hypothetical linker, which introduced its own complications. Nonetheless, Sr(NDC) provides clear evidence that the basic concept of the flexible slit pore can be engineered into a porous material. With further refinement and careful rational materials design, a similar material (with non-staggered slit pores and a larger probe occupiable pore volume than Sr(NDC) could provide deliverable capacity that approaches the ideal performance of the slit pore model and exceeds the best known benchmark materials. Concerted and extensive computational efforts will be needed to screen porous materials at the most accurate levels of theory and with fully flexible treatment to uncover them and motivate experimental testing.
Ultimately, this work outlines a less-explored path to discover the ultimate gas deliverable capacity materials via a concrete inverse design objective. Synthetic control over the chemical and structural properties of MOFs is always improving, and the structural motifs that give this adsorption model its outstanding deliverable capacities are increasingly being realized. For example, "adsorbaphores" (large conjugated linkers arranged in a parallel orientation) 70 and "nanographene" molecules 71 have been shown to provide advantageous adsorption properties when incorporated as MOF linkers. Further focus on these types of materials could lead to the discovery of a flexible slit pore geometry that provides the exact deliverable capacity advantages highlighted in this work. Most importantly, this adsorption model provides (1) a template for the rational materials design and (2) thermodynamic justification that the current methane deliverable capacity record holders could be surpassed by materials that undergo little to no volume change. This structural motif is not inherently limited to MOFs and should be explored/sought out in other systems known for gas storage capabilities, such as porous organic cages, covalent organic frameworks, etc.

Conflicts of interest
There are no conflicts to declare.
Materials Advanced Research Consortium (Hy-MARC) under Contract Number DE-AC04-94AL85000 with Sandia National Laboratories and Contract Number DE-AC02-05CH11231 with Lawrence Berkeley National Laboratory. Sandia National Laboratories is a multimission laboratory managed and operated by National Technology and Engineering Solutions of Sandia, LLC., a wholly owned subsidiary of Honeywell International, Inc., for the U.S. Department of Energy's National Nuclear Security Administration under contract DE-NA-0003525. We acknowledge the use of the Athena supercomputer through the HPC Midlands+ Consortium, and the ARCHER supercomputer through membership of the UK's HPC Materials Chemistry Consortium, which are funded by EPSRC Grants No. EP/P020232/1 and No. EP/R029431/1, respectively. We gratefully acknowledge Prof. Jeffrey R. Long for helpful discussions and measurements performed in his laboratory at Lawrence Berkeley National Laboratory.

Graphical TOC Entry
download file view on ChemRxiv FlexMethane.pdf (4.09 MiB) S1. Additional model details

S1.1. Effects of the E R F parameters on deliverable capacity
The shape of the framework potential, E F , impacts the adsorption energy landscape and thus the deliverable capacity of the adsorption model. However, significant variation in parameter selection for eq. (4) still leads to improved DC over the static slit pore. For example, a parameter set that produces an asymmetric single well potential still has favorable DC compared to the static slit pore. The parameters for these two variations are shown in Table S1, where θ ′ minimizes E R F and D i (θ ′ , L eq = 2σ) is the largest free sphere at that tilt angle. In other words, both materials have drastically improved DC over the static slit pore because their nonporous minimum energy states suppress adsorption at low pressures. The asymmetric single well potential is plotted in Figure S1c, and the asymmetric double well potential is plotted in Figure S1f. The isotherms resulting from these potentials are shown in blue in Figure S1. Notably, the asymmetric single well potential still demonstrates an inflection in the adsorption isotherm and maintains a larger DC than the optimized static slit pore (Figure 4).

S1.2. Cooperative adsorption effects from M > 1
The inflections in the isotherms in Figure 4 are dependent on a sufficiently large M . In other words, the porous state must have the capability to adsorb multiple methane molecules per slit pore in order to obtain the enhancement in deliverable capacity afforded by S-shaped isotherms. * mwitman@sandia.gov ; mdallen@sandia.gov  Figure S1: Isotherms plotted on a linear-linear scale (a,d), isotherms repeated on a log-log scale (b,e), and E R F (θ, L eq = 2σ) (c,f) for an asymmetric single well potential (a-c) and an asymmetric double well (d-f). Blue isotherms correspond to the he rotating slit pore isotherm, while dashed red lines correspond to the static slit pore with θ eq = 0 and θ eq = θ ′ .   Figure S2: Evolution of the isotherms and deliverable capacity as a function of M and L eq . The first column shows isotherms for various L eq , the second column shows the same data on a log-log scale, and the third column plots the deliverable capacity as a function of L eq . S3 S1.3. Adsorption volume discretization Figure S3 shows M = 4 discretization of channel sites within the slit pore. Each channel site has width, W = 2 1/6 σ G , and thus cross-sectional area, A site = W 2 . We denote the "free area", A o , as the area in which the particle can translate without overlapping a bouNDCry or an occupied cell. Overlap occurs when the particle's center is less than σ G from a site bouNDCry. This modeling of particle overlaps is important since such configurations are very high in energy and don't contribute to the partition function. This "free area" in lattice adsorption models is also commonly used to model dead volumes in porous materials. 1 In the case where there are no particles present (visualized by N = 0 in Figure S3), A o = (2W − σ G ) 2 /4. When 2 particles already occupy sites (blacked out cells for the N = 2 schematic in Figure S3), These free areas are schematically shown in green in Figure S3. Clearly A o depends on the number of particles adsorbed in the lattice model, but introduction of this complexity negates our ability to numerically integrate the partition function. Therefore we choose the intermediary A o value, which yields A o = 1.1 Å 2 when σ G = 3.73 Å. In other words, this method of calculating A o yields the highest fidelity at intermediate loadings, and subsequently overestimates the configurational entropy at high loadings and underestimates it a low loadings. Nonetheless, this discretization is a necessary step to make the model numerically integrable, to account for particle overlaps (i.e. much better than setting A o = A site ), and to allow the model results to be mapped to a volumetric basis for comparison with real materials.

S2.2. GCMC simulation details
GCMC simulations were performed using the RASPA code 14 with the TraPPE 15 force field Lennard-Jones parameters for methane and the UFF 16 force field Lennard-Jones parameters for framework atoms. Force field parameter values, structures, and MC sampling settings can be found in the supporting simulation files (S2.4).

S2.3. Isoreticular expansion of Sr(NDC)
More void volume can be imparted to the porous state of Sr(NDC) by performing isoreticular expansion with the hypothetical linkers (HYP, HYP2, HYP3) shown in Figure S4. Each ligand is of equal length (defined as distance between the carboxylic acid groups which form the 1-D SrO rods), but of varying width with W HYP < W HYP2 < W HYP3 . After constructing Sr(NDC) analogs using these hypothetical ligands via a 1-D rod MOF assembly algorithm, 17 we then applied the exact same DFT calculations as were described for Sr(NDC) in the manuscript. Visual comparison of the resulting F| DFT F and F| DFT F·G configurations are shown in Figure S5. The geometric differences between the F| DFT F and F| DFT F·G configurations between each linker significantly affect the adsorption properties. Again denoting θ ′ as the linker rotation angle in the F| DFT F state, Figure S5 shows that θ ′ HYP > θ ′ HYP2 > θ ′ HYP3 . This occurs because the wider linkers protrude further into the center of the bc face ( Figure S5), and steric repulsion with the opposite linker reduces how far each can rotate. A higher value of θ ′ reduces the porosity between the rotating pores, which helps boost DC. However, if the linker is not wide enough, a porous channel opens normal to the bc face, which also reduces the porosity of the F| DFT F configuration. This is especially noticeable in Sr(HYP). Therefore, a the optimal linker need have the perfect width to prevent adsorption between within the slit pore and prevent a channel from existing normal to the bc face.
None of the hypothetical linkers, however, achieve the perfect aspect ratio to reproduce the highly nonporous F| DFT F state of Sr(NDC). This can be seen from the adsorption isotherms performed on the F| DFT F state in Figure S6. With HYP, the large porous space in the bc face results in significant uptake in F| DFT F . HYP2 and HYP3 both show smaller, yet non-negligible uptake due to either too much void space in the middle of the bc face or too little linker rotation, respectively. Figure S6 also marks P = 35 bar, which is also sometimes selected as the P ads for technical storage targets. Despite the increase in porosity of the F| DFT F configurations, the maximum achievable deliverable capacity of the hypothetical structures with a P ads = 35 bar target could be competitive with the record-setting Co(BDP), which exhibits DC = 155 v(STP)/v. 18 Taking ρ| 35bar (F| DFT F·G ) − ρ| 5.8bar (F| DFT F ) to be the maximum possible deliverable capacity, Table S2 shows how the HYP2 and HYP3 ligands, in the best case scenario, could provide comparable performance to Co(BDP). If a structure could be optimized that, like Sr(NDC), has a more nonporous F| DFT F state, it would then have the potential to exceed Co(BDP) performance with a P ads = 35 bar target.

S2.4. Supplementary files
The following structural and simulation files are provided in the supplementary_files.zip, as well as the IPython notebook for the slit pore model:

S3. Sr(NDC) synthesis and characterization
Sr(NDC) was synthesized via a solvothermal reaction by reacting strontium nitrate with 1,4-naphtalenenaphthalenedicarboxylic acid (1,4-H 2 NDA) using a procedure previously described by Raja et al. 13 In order to obtain sufficient powder for hydrogen isotherm measurements, the synthesis was scaled up by a factor of 50 compared to the literature method. This was achieved by conducting the solvothermal synthesis in five 250 mL glass bottles and combining the product. In each 250 mL bottle, 0.43 g of 1,4-H 2 NDA was dissolved in 70.0 mL DMF, 20.0 mL EtOH, and 10.0 mL of H 2 O. To this solution, 0.85 g of Sr(NO 3 ) 2 was added and the reaction mixture was stirred at room temperature for 30 min until all the solids completely dissolved. The solution was then placed in an oven at 120 • C for 48 hours. The resulting rod-shaped crystals were washed using 100 mL of EtOH and dried overnight under ambient conditions, then dried in an oven at 80 • C for 2 hours. Total final yield of the dried product was 4.524 g. The activation of Sr(NDC) was conducted prior to hydrogen isotherm measurements by heating the sample for 24 hours at 300 • C in vacuum. The materials were characterized using in situ XRD on a Oxford Diffraction Supernova diffractometer in capillary mode using Cu Kα radiation. The samples were loaded into 0.7 mm diameter capillaries inside a glovebox. Measurements were recorded using a CCD detector placed at 77 mm from the samples with an exposure time of 60 seconds. The recorded 2D images were added and integrated to generate a 1D pattern.