Adsorption and dynamics of linear and mono-branched hexane isomers in MIL-140 metal–organic frameworks

Abstract

Recent breakthrough experiments revealed the iso-reticular Zr-MOFs, MIL-140B and MIL-140C, as promising sorbents for the separation of C6 isomers. Interestingly while the ultra-small pore MIL-140B exhibited hexane isomer sorption hierarchy according to the normal boiling point order (n-C6 > 3MP (3-methyl pentane)), an uncommon shift in the elution order was observed in the larger pore MIL-140C. It was only speculated that the flexibility of the MOFs might be the origin of this intriguing behavior. Herein, flexible force field hybrid osmotic Monte Carlo combined with molecular dynamics simulations were carried out to unravel the microscopic mechanism of the adsorption and dynamics of both C6 isomers in MIL140B and MIL140C. Thermodynamically preferred adsorption of n-C6 over 3MP was predicted for MIL-140B and to a slightly less extent for MIL-140C. Interestingly while the mobility of n-C6 was found to remain higher than that of 3MP in the whole range of loading for MIL-140B, 3MP becomes more mobile than n-C6 at saturation in MIL-140C. This suggests that this kinetics order is most probably the origin of the inversion of the elution order observed experimentally for MIL-140C. The translational and rotational dynamics of the two guests in MIL-140B and MIL-140C was further understood in-depth.

---

1 Introduction

The separation of hexane isomers for the enrichment of the octane number (RON) of gasoline is a key process in the petroleum industry.^1,2 A family of metal organic frameworks (MOFs)^3–9 including Fe₂(BDP)₃,¹⁰ Zr-abtc,¹¹ MIL-140B,¹² MIL-53(Fe)-(CF₃)₂,¹³ Ca(H2tcpb),¹⁴ Al-bttotb¹⁵ and more recently MIL-160(Al)¹⁶ has been revealed to be alternative sorbents to the standard zeolites^17–19 for discriminating the hexane isomers into valuable fractions according to the degree of branching via thermodynamically controlled, kinetically driven or molecular sieving-controlled mechanisms. In particular, the pore dimension of an isoreticular series of channels like Zr-MOFs MIL-140A (3.2 Å), MIL-140B (4.0 Å) and MIL-140C (5.7 Å)²⁰ (Fig. 1a and b) has been shown to be highly modulable by changing the nature/length of the linkers (R = C₆H₄, C₁₀H₆, and C₁₂H₈ for MIL-140A, MIL-140B and MIL-140C, respectively) enabling a fine tuning of their hexane isomer separation properties.^12,21 Indeed, MIL-140B and MIL-140C were identified as the most attractive sorbents for the selective separation of hexane isomers while MIL-140A was not appropriate because its pore size was too small (3.2 Å) to accommodate most of the C6 molecules. It was demonstrated that MIL-140B can separate the hexane isomers in the order of linear (n-C6) > mono-branched (3-methylpentane (3MP)) > di-branched at 343 K together with a high selectivity up to 10 as well as a sorption hierarchy similar to the normal boiling point order of the compounds.²² In contrast, a shift in the elution order occurring between n-C6 and 3MP was observed in the experimental breakthrough curve obtained for MIL-140C.¹² This work aims to unravel the microscopic mechanism at the origin of this uncommon reverse elution order n-C6 and 3MP between MIL-140C and MIL-140B using flexible force field hybrid osmotic Monte Carlo and molecular dynamics simulations.

Fig. 1 Illustration of the channel-like structure of MIL-140B (a) and MIL-140C (b) MOF materials. Gray, red, white and pink colors correspond to carbon, oxygen, hydrogen and aluminium atoms, respectively.

2 Computational details

Single component and binary mixture adsorption isotherms for n-C6 and 3MP isomers in MIL-140B and MIL-140C were predicted by using hybrid osmotic Monte Carlo (HOMC) simulations. The HOMC simulations are based on the introduction of a molecular dynamic (MD) step in a Monte Carlo (MC) scheme through the osmotic statistical ensemble.²³ Theoretical background ruling the HOMC method is detailed in ref. 23 and 24. In the HOMC simulation step, the MD part consisted of a run of 10⁵ steps using a time-step of 0.001 ps i.e. 100 ps in order to make sure that the thermodynamic equilibrium was achieved.²⁴ In the grand canonical Monte Carlo (GCMC) simulations, the chemical potential, temperature, and volume of the system were fixed. All simulations consisted of 10⁵ MC steps. Insertion/deletion, translation, rotation and internal rotation trial moves have been considered with the following frequency: 0.4, 0.2, 0.2 and 0.2. For the gas mixture simulations, additional swapping trial move was considered with a frequency of 0.1 and of 0.3 for the insertion/deletion move. The MD steps were performed in the NσT (N: number of particles, σ: constraint and T: temperature = 343 K) ensemble where the thermostat (at a relaxation time of 0.1 ps) and the barostat (at a relaxation time of 0.5 ps) were both considered by means of the Nose–Hoover algorithms.^25,26 Each MD step was accepted following the metropolis criterion.²³

These simulations used a simulation box corresponding to a 1 × 2 × 4 supercell of the MIL-140B and MIL-140C crystal structure reported previously.²⁰ The C6 isomers, n-C6 and 3MP, were modeled by means of the flexible OPLS all atoms force field.²⁷ To be in line with the previous studies on MIL-140B,²¹ the UFF²⁸ (for the inorganic nodes) and DREIDING²⁹ (for the organic linkers) force fields were combined to describe the flexibility of the MIL-140 materials. This combination of the two force fields treats the MIL-140 flexibility with bonding, bending, torsion and improper intra-molecular potential terms while non-bonded interactions are described with Lennard-Jones (LJ) and electrostatic potentials. Preliminary molecular dynamics (MD) simulations on MIL-140C led to simulated unit cell parameters (a = 30.62 Å, b = 15.49 Å and c = 7.31 Å, β = 91.02°) in fair concordance with the corresponding experimental data²⁰a = 31.03 Å, b = 15.51 Å, c = 7.82 Å, and β = 93.26°. Pore size distribution for both MIL-140B and MIL-140C was calculated from the Zeo++ software³⁰ based on the Voronoi decomposition, which provides a graphical representation of the void space for a given arrangement of atoms in a periodic domain. The resulting Voronoi network was analyzed to obtain the diameter of the largest included sphere frequently used to describe pore geometry. The simulated PSD profiles were found to be centered around the mean value reported previously for the crystal structure of the two MOFs. This favorable comparison highlights that the LJ parameters as well as the intra-molecular terms are reliable to describe the MIL-140 frameworks. The cross C6/MIL-140 LJ interaction parameters were calculated using the Lorentz-Berthelot mixing rule.^31,32 All calculated LJ interactions were truncated with the use of a cutoff distance of 12 Å. Electrostatic interactions were modeled by means of the Ewald summation technique³³ using a convergence parameter of 10⁻⁶.

In complement to these HOMC simulations, additional MD simulations were carried out in the NσT ensemble where σ is the constraint (1 bar) to calculate the mean square displacement (MSD) of the hexane isomers to explore the sorption hierarchy from a dynamic point of view.

with r_com,i is the position of the centre of mass of molecule i, t₀ is the time origin, N is the number of molecules and N₀ is the number of t₀.

3 Breakthrough experiments

To collect the single component adsorption equilibrium isotherms of 3MP and n-C6, breakthrough experiments were performed in a stainless steel fixed bed column (internal diameter 4.6 mm and length 100 mm) containing 590 mg of MIL-140. The MIL-140B and MIL-140C samples were taken from the batch synthesized and characterized previously.¹² The feed is set up using a mass flow controller (helium) and a syringe pump (paraffins), with the mixture being completely vaporized in a pre-heater before entering the column. The outlet stream of the fixed bed is directed to a thermal conductivity detector (TCD), which continuously measures the concentration profiles of the hydrocarbon species as a function of time. The equilibrium loading (amount adsorbed) of the components is obtained by integrating the concentration profiles (molar rate history) of the breakthrough curves.

4 Results and discussion

The single component adsorption isotherms calculated for n-C6 and 3MP in MIL-140B and MIL-140C at 343 K are reported in Fig. 2a. MIL-140B is shown to have a higher affinity for n-C6 than 3MP at the initial stage of adsorption (steeper adsorption isotherm) consistent with a higher adsorption enthalpy calculated at low loading for (n-C6@MIL-140B) = −84.5 kJ mol⁻¹vs. >Q_st(n-C6@MIL-140C) = −71.5 kJ mol⁻¹ (Q_st = RT − [〈UN〉 − 〈U〉〈N〉]/[〈N²〉 − 〈N〉²] where U, N, R, and T are the host–guest interactions, the number of molecules, the perfect gas constant and the temperature, respectively). We equally predict a larger n-C6 amount adsorbed in line with the linear shape of this isomer that accommodates better the shape of the small triangular channels of the MOF (about 4 Å as depicted by the pore size distribution (PSD) of both empty MOFs calculated from MD simulations for the empty structures, Fig. 2b) as compared to the bulkier 3MP. In the case of MIL-140C, the amount of n-C6 and 3MP adsorbed is relatively similar with an affinity for n-C6 still slightly higher as also revealed by a higher adsorption enthalpy at low coverage calculated for the linear isomer (−71.5 kJ mol⁻¹vs. −68.7 kJ mol⁻¹). Indeed, the larger pore size of MIL-140C (about 5.7 Å as illustrated in Fig. 2b enables it to host more linear and branched alkanes in the MOF channel leading to adsorption uptakes larger than in MIL-140B. These single component adsorption isotherms simulated for both MIL140B and MIL140C are in good agreement with the corresponding experimental data shown in Fig. 2a.

Fig. 2 (a) HOMC calculated single component adsorption isotherms of hexane isomers at 343 K for (a) MIL-140B and MIL-140C MOF materials. n-C6 in MIL-140B (), 3MP in MIL-140B (), n-C6 in MIL-140C () and 3MP in MIL-140C () and experimental data n-C6 in MIL-140B (), 3MP in MIL-140B (), n-C6 in MIL-140C () and 3MP in MIL-140C (). (b) Pore size distribution of both empty MIL-140B (filled black histogram) and MIL-140C (empty red histogram) MOF materials. HOMC calculated C6/3MP binary adsorption isotherms at 343 K for (c) MIL-140B and (d) MIL-140C MOF materials, n-C6 (), 3MP () and total ().

The n-C6/3MP binary mixture adsorption isotherms for both MIL-140B and MIL-140C materials are shown in Fig. 2c and d, respectively. These simulated data highlight that these two MOFs show the same n-C6 > 3MP adsorption hierarchy following the same trend as in the single component adsorption isotherms. Let us mention that this trend is also observed at low vapour pressure (Fig. S1 of ESI†). The preferential adsorption of n-C6 over 3MP is more pronounced at low loading for MIL-140B in line with a more pronounced n-C6/3MP adsorption enthalpy difference discussed above. Interestingly, the predicted thermodynamic sorption hierarchy for MIL-140B is in line with the experimental hierarchy determined previously by breakthrough curve measurements performed on the quinary C6, 3MP, 2MP, 2,2-DMB and 2,3-DMB (DMB = DiMethyl-Butane) mixtures with a higher sorption capacity along with a longer retention time for n-C6 vs. 3MP. In contrast, the opposite trend is observed for MIL-140C, the same breakthrough curve experiments revealing that 3MP is the last component to saturate in the column and indeed there is a shift between 3MP and n-C6 in the order of elution, the sorption hierarchy being 3MP > n-C6.

This deviation between thermodynamic adsorption isotherms and breakthrough experiments may suggest that there is a kinetics contribution to this dynamic separation in these two channel like MOFs characterized by a certain degree of flexibility of the organic linkers. As a further step, MD simulations were carried out at low (1 kPa) and high (25 kPa) vapour pressures to explore the dynamics of the two guest molecules as single components in MIL-140B and MIL-140C materials. The corresponding mean square displacements (MSD) of the two isomers averaged over multiple time origins and several MD trajectories for both MOFs are reported in Fig. 3. As shown in Fig. 3a, the mobility sequence of the hexane isomers at very low loading in the MIL-140B is n-C6 > 3MP. Since 3MP is a bulkier molecule than n-C6, its dynamic is slowed down in the narrow channel of MIL-140B. In the case of MIL-140C, its larger pore dimension enables both isomers to adapt a similar kinetics. At high vapor pressure, the dynamic trend remains the same in MIL-140B, the calculated MSD for n-C6 being still higher than the value obtained for 3MP. Indeed there is no kinetics limitation that would hamper the thermodynamics sequence. This trend clearly highlights that the higher affinity predicted for n-C6 vs. 3MP governs the higher sorption elution for the linear alkane as well as the longer retention time, the host/guest interactions governing the breakthrough sequence.

Fig. 3 Mean square displacement of single components n-C6 and 3MP in both MIL-140B and MIL-140C materials at 343 K for a vapour pressure (a) 0.1 kPa and (b) 25 kPa, n-C6 in MIL-140B (—), 3MP in MIL-140B (––), n-C6 in MIL-140C () and 3MP in MIL-140C ().

Regarding MIL-140C, while the mobilities of n-C6 ∼ 3MP are similar at low loading (Fig. 3a), interestingly 3MP becomes more mobile than n-C6 at higher loading (Fig. 3b). Indeed this pronounced decrease of the mobility of n-C6 with respect to 3MP probably contributes to counterbalance the host/guest interactions in favor of C6 that could explain that both hexane isomers elute practically at the same time and that the sorption hierarchy is in favor of 3MP.

To further gain microscopic insight into the behavior of the two isomers in both MOFs, two dimensional (2D) densities of the center of mass of each single component C6 molecule according to the yz plan with z is the channel direction (Fig. 1) were calculated at saturation conditions. Fig. 4 reports the 2D density of n-C6 confined in both MIL-140B and MIL-140C along one representative channel. As shown in Fig. 4a, the guest population is well ordered in the MIL-140B channel with spatial patterns spaced of 2 Å suggesting a single file diffusion process with a short translational jump. Single file diffusion was corroborated by the calculation of coefficient α in the generalized Einstein's equation;^34,35 MSD = At^α where t is the time and A is a coefficient. For α = 1, α < 1 and α > 1 the dynamics is diffusive, sub-diffusive or super-diffusive, respectively. For α = 1/2, the dynamics corresponds to the single file diffusion.³⁴ In the case of n-C6 confined in MIL-140B α = 0.59 close to 0.5 that corroborates a single file diffusion. At low loading (0.1 kPa) (Fig. 3), the probability of finding two molecules in the same channel is very low. We typically have one molecule per channel, which explains the absence of single-file diffusion. At this filling rate, we have a diffusive regime instead. That was corroborated from the calculation of the alpha coefficient in the generalized Einstein relation close to 1.0 (0.95) in line with a diffusive dynamic.

Fig. 4 Two dimensional density of centres of mass of (a) n-C6 in molecules MIL-140B, (b) 3MP molecules in MIL-140B, (c) n-C6 molecules in MIL-140C and (d) 3MP molecules in MIL-140C. These data are analyzed from the MD simulations performed at 343 K and a vapour pressure of 25 kPa.

Fig. 4b exhibits a fixed position of 3MP molecules into the channel of MIL-140B highlighting small translational degrees of freedom in line with the slow dynamics reported in Fig. 3b. In contrast, Fig. 4c reveals that in MIL-140C the spatial patterns are separated by a distance about 7 Å too large to be considered as a translational jump that suggests fixed positions of n-C6 molecules into the channel with few translational degrees of freedom that is in good agreement with a slow single file dynamics of this molecule reported in Fig. 3b. n-C6 optimizes its spatial organization into the channel leading to a decrease in conformational entropy. Fig. 4d shows that the 3MP molecules are not less anchored in the channel of MIL-140C highlighting an increase in the translational degree of freedom in comparison with n-C6 which is in good agreement with the reported MSD in Fig. 3b.

To complete the investigation of the translational motions, rotational degrees of freedom were examined by evaluating the distribution of the angle (Θ_z) between the unit vector along the z direction and the vector related to the end-to-end distance of C6 isomer molecules (d_ee). As shown in Fig. 5a, n-C6 and 3MP molecules present an elongated structure in the MIL-140B material because Θ_z = 0° and 180°. This is the result of the high confinement effect where only elongated conformation is allowed. This point is corroborated by the calculation of the gyration radius according to the three directions (R_gx, R_gy, R_gz). Gyration radii were calculated from the inertia tensor and its diagonalisation.³⁶ For n-C6 in MIL-140B, R_gz = 4.3 Å against 2.7 Å in MIL-140C that highlights a more elongated shape in MIL-140B. In MIL-140C, the angular situation differs for both C6 isomers. Indeed, 3MP molecules adopt angles between 140° and 180° and between 0° and 40° due to the increase in the pore size of MIL-140C and the flexibility of 3MP. Indeed, Fig. 5b which reports the dihedral angles distribution of C6 isomers shows that 3MP is more flexible than n-C6 in both MIL-140B and MIL-140C. Conversely, n-C6 shows a bimodal distribution of Θ_z centered on 40° and 180° highlighting two preferential conformations. These results emphasize that n-C6 is strongly anchored in MIL-140B with two possible conformations.

Fig. 5 Distribution of (a) angle defined at the top of the figure between the unit vector along z direction and the vector related to the end-to-end distance of hexane isomer molecules (d_ee) and (b) the dihedral angles in hexane isomer molecules as illustrated at the top of the figure at 343 K and for a vapour pressure of 25 kPa, n-C6 in MIL-140B (—), 3MP in MIL-140B (––), n-C6 in MIL-140C () and 3MP in MIL-140C ().

As a further stage, the z-positions of three tagged n-C6 molecules confined in one channel of MIL-140C were followed as a function of time. As depicted in Fig. 6a, diffusion of 3.5 Å is observed around 4 ns for 3 molecules in a concerted dynamics. Furthermore, as shown in Fig. 6b, these translational diffusions are accompanied by a concerted rotational jump of 100° from 40° to 140° and from 140° to 40°. The dynamics of n-C6 confined in MIL-140B is then controlled by collective and concerted translational and rotational movements. MSD and Θ_z distributions of C6 and 3MP in the binary mixtures were calculated at the saturation to examine the influence of 3MP molecules on the anchoring of n-C6 in the case of pure components. Fig. 7a shows that n-C6 molecules adopt any specific orientations in a similar way that in the single component (Fig. 5a) contrary to the 3MP where all orientations are sampled in both binary mixtures and pure components. This result suggests that 3MP molecules prevent n-C6 packing and tend to break a strict ordering of the molecules along the channel. This is corroborated by the calculation of the two dimensional density distributions of two n-C6 molecules as in Fig. 7b and c where an absence of preferential locations is highlighted. These results allow us to well apprehend the observed difference between single and mixture adsorption isotherms that are more significant in the case of MIL-140C since a decrease of 86% is observed. That seems to be the result of favorable packing of n-C6 within MIL-140C minimizing the adsorption of 3MP that could harm n-C6 packing and break a strict ordering of the molecules along the channel. Interestingly, this entropic consideration was also mentioned to explain the difference between single and mixture isotherms of isomer alkane adsorption in silicalite.^37,38

Fig. 6 Time evolution of (a) z-position of three n-C6 molecules confined in one channel, molecule 1 (—), molecule 2 (––), molecule 3 () and (b) angle between the unit vector along y direction and the vector related to the end-to-end distance of n-C6 molecules in MIL-140B at 343 K and for vapour pressure of 25 kPa. On the right of (b), snapshots of the molecules are provided between 0 and 4 ns and between 4 and 10 ns.

Fig. 7 (a) Distribution of the angle (Θ_z) of n-C6 (—) and 3MP (––) in the binary mixture at 343 K and for a total vapour pressure of 25 kPa. (b and c) Two dimensional density of centres of mass of two different n-C6 molecules in MIL-140C.

5 Conclusion

In this work, hybrid osmotic Monte Carlo and Molecular Dynamics simulations have been carried out to reveal the molecular origin of the separation of linear to mono-branched hexane isomers in the two Zr-MOFs, MIL-140B and MIL-140C. Single component adsorption isotherms were first simulated, and they were found to be in fair agreement with the corresponding experimental data. Regarding the binary mixture, we found a thermodynamic sorption hierarchy n-C6 > 3MP in both MIL-140B and MIL-140C materials with a selectivity between 5 and 7, respectively, due to the most favorable MIL-140/n-C6 interactions. Although this sequence is in line with the elution order obtained experimentally from the breakthrough curve for MIL-140B, there is an inversion for MIL-140C (3MP > n-C6) suggesting that the separation is not solely controlled by the thermodynamics. The dynamics of the single component C6 isomers was further investigated using MD simulations at saturation that revealed a substantial slowing down of the mobility of n-C6 confined in MIL-140C compared to 3MP that was ascribed to strong packing of n-C6 molecules. For confined n-C6 molecules in MIL-140C, individual translational jumps have been evidenced. These simulations helped us to better understand the different breakthrough curve trends observed for MIL-140B and MIL-140C. In MIL-140B, the separation is likely to be predominantly controlled via thermodynamics with the host–guest interactions in favor of n-C6 leading to an elution order of n-C6 > 3MP. In MIL-140C, n-C6 has been found to be still slightly more preferentially adsorbed than 3MP from a thermodynamics point of view however compensated by a dynamics of n-C6 that becomes slower than that of n-C6 resulting into an elution of n-C6 and 3MP practically at the same time.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

G. M. thanks ANR MEACOPA (ANR-17-CE29-0003) for funding. J. A. C. S. thanks the Foundation for Science and Technology (FCT, Portugal) for financial support through national funds FCT/MCTES (PIDDAC) to CIMO (UIDB/00690/2020 and UIDP/00690/2020) and SusTEC (LA/P/0007/2020).

Notes and references

1. N. Cusher, RON numbers in handbook of Petroleum Refinning Processes, ed. R. A. Meyers, 3rd edn, 2004.
2. R. A. Myers, Handbook of Petroleum Refining Processesed, McGraw-Hill, NewYork, 2004.
3. B. Chen, C. Liang, J. Yang, D. Contreras, Y. Clanct, E. Lobkovsky, O. Yaghi and S. Dai, Angew. Chem., Int. Ed., 2006, 45, 1390.
4. P. Barcia, F. Zapata, J. Silva, A. Rodrigues and B. Chen, J. Phys. Chem. B, 2007, 111, 6101.
5. D. Dubbeldam, C. Galvin, K. Walton, D. Ellis and R. Snurr, J. Am. Chem. Soc., 2008, 130, 10884.
6. P. Barcia, D. Guimaraes, P. Mendes, J. Silva, V. Guillerm, H. Chevreau, C. Serre and A. Rodrigues, Microporous Mesoporous Mater., 2011, 139, 67.
7. Y. Ling, Z.-X. Chen, F.-P. Zhai, Y.-M. Zhou, L.-H. Weng and D.-Y. Zhao, Chem. Commun., 2011, 47, 7197.
8. D. Dubbeldam, R. Krishna, S. Calero and A. O. Yazaydin, Angew. Chem., Int. Ed., 2012, 51, 11867.
9. T. Duerinck and J. Denayer, Adsorption, 2014, 20, 251.
10. Z. Herm, B. Wiers, J. Mason, J. van Baten, M. Hudson, P. Zajdel, C. Brown, N. Masciocchi, R. Krishna and J. Long, Science, 2013, 340, 960–964.
11. H. Wang, X. Dong, J. Lin, S. J. Teat, S. Jensen, J. Cure, E. V. Alexandrov, Q. Xia, K. Tan, Q. Wang, D. H. Olson, D. M. Proserpio, Y. J. Chabal, T. Thonhauser, J. Sun, Y. Han and J. Li, Nat. Commun., 2018, 9, 1745.
12. A. Henrique, T. Maity, H. Zhao, P. Brantuas, A. Rodrigues, F. Nouar, A. Ghoufi, G. Maurin, J. Silva and C. Serre, J. Mater. A, 2020, 8, 17780.
13. P. A. P. Mendes, P. Horcajada, S. Rives, H. Ren, A. E. Rodrigues, T. Devic, E. Magnier, P. Trens, H. Jobic, J. Ollivier, G. Maurin, C. Serre and J. A. C. Silva, Adv. Funct. Mater., 2014, 24, 7666–7673.
14. H. Wang, X. Dong, E. Velasco, D. Olson, Y. Han and J. Li, Energy Environ. Sci., 2018, 11, 1226–1231.
15. L. Yu, X. Dong, Q. Gong, S. Acharya, Y. Lin, H. Wang, Y. Han, T. Thonhauser and J. Li, J. Am. Chem. Soc., 2020, 142, 6925–6929.
16. P. F. Brantuas, A. Henrique, M. Wahiduzzaman, A. vib Wedelstedt, T. Maity, A. E. Rodrigues, F. Nouar, U.-H. Lee, K. H. Cho, G. Maurin, J. A. C. Silva and C. Serre, Adv. Sci., 2022, 9, 2201494.
17. E. Jolimaitre, K. Ragil, M. Tayakout-Fayolle and C. Jallut, AIChE J., 2002, 48, 1927–1937.
18. R. Krishna and J. van Baten, Sep. Sci. Technol., 2007, 55, 246–255.
19. P. Barcia, J. Silva and A. Rodriguez, AIChE J., 2007, 53, 1970–1981.
20. V. Guillerm, F. Ragon, M. Dan-Hardi, T. Devic, M. Vishnuvarthan, B. Campo, A. Vimont, G. Clet, Q. Yang, G. Maurin, G. Férey, A. Vittadini, S. Gross and C. Serre, Angew. Chem., Int. Ed., 2012, 51, 9267.
21. H. Zhao, G. Maurin and A. Ghoufi, J. Phys. Chem. C, 2022, 126, 2905–2911.
22. G. Cholakov, W. Wakeham and R. Stateva, Fluid Phase Equilib., 1999, 163, 21–42.
23. A. Ghoufi and G. Maurin, J. Phys. Chem. C, 2010, 114, 6496–6502.
24. H. Zhao, G. Maurin and A. Ghoufi, J. Chem. Phys., 2021, 154, 084702.
25. S. Nose, J. Chem. Phys., 1984, 81, 511.
26. W. Hoover, Phys. Rev. A, 1985, 31, 1695.
27. W. Jorgensen, D. Maxwell and J. Tirado-Rives, J. Am. Chem. Soc., 1996, 118, 11225.
28. A. Rappe, C. Casewit, K. Colwell, W. Goddard III and W. Skiff, J. Am. Chem. Soc., 1992, 114, 10024–10035.
29. S. Mayo, B. Olafson and W. Goddard, J. Phys. Chem., 1990, 94, 8897–8909.
30. T. Willems, C. Rycroft, M. Kazi, J. Meza and M. Haranczyk, Microporous Mesoporous Mater., 2012, 149, 134–141.
31. H. A. Lorentz, Ann. Phys., 1881, 248, 127.
32. D. Berthelot, C. R. Hebd. Seances Acad. Sci., 1898, 1703.
33. P. Ewald, Ann. Phys., 1921, 64, 253–287.
34. A. Ghoufi and G. Maurin, J. Phys. Chem. C, 2019, 123, 17360–17367.
35. R. Boulé, C. Roiland, T. Bataille, L. L. Pollès, N. Audebrand and A. Ghoufi, J. Phys. Chem. Lett., 2019, 10, 1698–1708.
36. A. Mendes, C. Bonal, N. Morel-Desrosier, J. Morel and P. Malfreyt, J. Phys. Chem. B, 2002, 106, 4516.
37. T. Vlugt, R. Krishna and B. Smit, J. Phys. Chem. B, 1999, 103, 1102–1118.
38. T. Titze, C. Chmelik, J. Karger, J. van Baten and R. Krishna, J. Phys. Chem. C, 2014, 118, 2660–2665.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: https://doi.org/10.1039/d2cp05371c

---