Membrane separation study for methane-hydrogen gas mixtures by molecular simulations

TL;DR

This study uses molecular dynamics simulations to analyze methane-hydrogen separation in silica zeolite membranes, revealing that permeation selectivities favor methane and differ from equilibrium adsorption predictions.

Contribution

It introduces a hybrid non-equilibrium molecular dynamics methodology to simulate gas transport in real zeolite structures, providing detailed insights into permeation selectivities.

Findings

Permeation selectivities favor methane over hydrogen.

Transport selectivities differ from equilibrium adsorption ratios.

Selectivity depends on membrane loading and component mobilities.

Abstract

Direct simulation results for stationary gas transport through pure silica zeolite membranes (MFI, LTA and DDR types) are presented using a hybrid, non-equilibrium molecular dynamics simulation methodology introduced recently. The intermolecular potential models for the investigated CH$_{4}$ and H$_{2}$ gases were taken from literature. For different zeolites, the same atomic (Si and O) interaction parameters were used, and the membranes were constructed according to their real (MFI, LTA, or DDR) crystal structures. A realistic nature of the applied potential parameters was tested by performing equilibrium adsorption simulations and by comparing the calculated results with the data of experimental adsorption isotherms. The results of transport simulations carried out at 25$^0$C and 125$^0$C, and at 2.5, 5 or 10 bar clearly show that the permeation selectivities of CH$_{4}$ are higher than the corresponding permeability ratios of pure components, and significantly differ from the equilibrium selectivities in mixture adsorptions. We experienced a transport selectivity in favor of CH$_{4}$ in only one case. A large discrepancy between different types of selectivity data can be attributed to dissimilar mobilities of the components in a membrane, their dependence on the loading of a membrane, and the unlike adsorption preferences of the gas molecules.

Tables

Table 1. Table 1: Thicknesses and surface sizes of the membranes used in the transport simulations.

		zig-zag (or rotated) channels	straight channels
MFI	thickness/nm	2.4200	2.2300
	surface size/nm2	41.5035	42.2436
LTA	thickness/nm	2.6440	2.6351
	surface size/nm2	49.6903	50.2267
DDR	thickness/nm	2.6494	3.0200
	surface size/nm2	44.9719	48.6832

Table 2. Table 2: Effective pair potential parameters (size and energy) used in the simulations.

Atom/molecule	$\sigma$ / nm	$(ϵ/k_B$) / K
O [23]	0.270	130.0
Si [25]	0.070	20.0
CH4 [23]	0.373	147.9
H2 [24]	0.296	36.7
Soft repulsive wall	0.300	120.0

Table 3. Table 3: Equilibrium and dynamic concentrations ( c 𝑐 c ) on the different zeolites [adsorbent (denoted by a) or membranes (denoted by m)] and transport fluxes ( j 𝑗 j ) obtained with equimolar CH 4 -H 2 gas mixture.

	zig-zag (or rotated) channels				straight channels
$T$ / K	298.15		398.15		298.15		398.15
$p$ / bar	2.5	10	2.5	10	2.5	10	2.5	10
$c_{{\mathrm{CH}}_4}$/(mol/dm3 MFIa)	1.255	2.887	0.224	0.787	1.255	2.887	0.224	0.787
$c_{H_2}$/(mol/dm3 MFIa)	0.027	0.065	0.016	0.058	0.027	0.065	0.016	0.058
$c_{{\mathrm{CH}}_4}$/(mol/dm3 MFIm)	0.848	1.493	0.256	0.455	1.007	1.747	0.324	0.524
$c_{H_2}$/(mol/dm3 MFIm)	0.020	0.054	0.010	0.034	0.018	0.048	0.008	0.030
$j_{{\mathrm{CH}}_4}$/(ns$\cdot$nm2 MFIm)	0.17	1.00	0.08	0.37	0.09	0.94	0.14	0.66
$j_{H_2}$/(ns nm2 MFIm)	0.16	0.51	0.09	0.51	0.10	0.78	0.17	1.15
$c_{{\mathrm{CH}}_4}$/(mol/dm3 LTAa)∗0.6	0.47	1.19	0.10	0.34	0.47	1.19	0.10	0.34
$c_{H_2}$/(mol/dm3 LTAa)	0.039	0.131	0.023	0.087	0.039	0.131	0.023	0.087
$c_{{\mathrm{CH}}_4}$/(mol/dm3 LTAm)	0.441	0.918	0.140	0.286	0.637	1.103	0.196	0.358
$c_{H_2}$/(mol/dm3 LTAm)	0.024	0.075	0.012	0.045	0.019	0.064	0.010	0.036
$j_{{\mathrm{CH}}_4}$/(ns$\cdot$nm2 LTAm)	0.14	0.56	0.08	0.36	0.17	0.62	0.10	0.37
$j_{H_2}$/(ns$\cdot$nm2 LTAm)	0.15	0.90	0.11	0.80	0.15	0.92	0.12	0.74
$c_{{\mathrm{CH}}_4}$/(mol/dm3 DDRa)∗0.6	1.58	2.43	0.31	0.93	1.58	2.43	0.31	0.93
$c_{H_2}$/(mol/dm3 DDRa)∗0.7	0.02	0.05	0.01	0.05	0.02	0.05	0.01	0.05
$c_{{\mathrm{CH}}_4}$/(mol/dm3 DDRm)	0.808	1.484	0.294	0.435	0.953	1.486	0.296	0.473
$c_{H_2}$/(mol/dm3 DDRm)	0.016	0.048	0.007	0.025	0.016	0.048	0.007	0.025
$j_{{\mathrm{CH}}_4}$/(ns$\cdot$nm2 DDRm)	0.06	0.11	0.01	0.08	0.03	0.10	0.02	0.06
$j_{H_2}$/(ns$\cdot$nm2 DDRm)	0.06	0.34	0.04	0.35	0.06	0.35	0.03	0.30
∗ due to inaccessible cavities, simulated data are corrected by the factor indicated.

Equations

$$\left|\frac{N_{\mathrm{control}\mathrm{cell}}\pm 1}{N_{\mathrm{control}\mathrm{cell}}}{p}_{\mathrm{control}\mathrm{cell}}-p_{\mathrm{target}}\right|<\left|p_{\mathrm{control}\mathrm{cell}}-p_{\mathrm{target}}\right|,$$

$$u^{\alpha\beta}(r)=\left\{\begin{array}[]{ll}u_{\mathrm{LJ}}^{\alpha\beta}(r)-u_{\mathrm{c}}^{\alpha\beta}&\;\mbox{for}\;\;r<r_{\mathrm{c}}^{\alpha\beta},\\ 0&\;\mbox{for}\;\;r\geqslant r_{\mathrm{c}}^{\alpha\beta},\end{array}\right.$$

$$u_{\mathrm{LJ}}^{\alpha \beta }(r)=4{ϵ}^{\alpha \beta }\left[{\left(\frac{{\sigma }^{\alpha \beta }}{r}\right)}^{12}-{\left(\frac{{\sigma }^{\alpha \beta }}{r}\right)}^6\right]$$

$$S_{\mathrm{E}}=\frac{q_{{\mathrm{CH}}_4}}{q_{{\mathrm{H}}_2}}.$$

$$S_{\mathrm{P}}=\frac{j_{{\mathrm{CH}}_4}}{j_{{\mathrm{H}}_2}},$$

Full text

Key words: gas permeation, zeolite membrane, steady-state, molecular dynamics

PACS: 02.70.Ns, 07.05.Tp, 68.43.-h, 68.43.Jk

Ключовi слова: ïðîíèêíåííÿ ãàçó, öåîëiòîâà ìåìáðàíà, ñòàöiîíàðíèé ðåæèì, ìîëåêóëÿðíà äèíàìiêà

1 Introduction

Zeolites are made up of silicon, aluminum, and oxygen atoms linked together so that they form structurally well-defined pores. A high regularity of the structure distinguishes zeolites from other porous materials and makes their high selectivity in catalytic and separation processes possible [1, 2]. A usual size of zeolite micropores is similar to that of many small molecules. Therefore, it is possible for some molecules to enter the zeolite and then get stuck in the pores, where they can even react with each other, while the other molecules can move through the zeolite channels faster. In separation processes, the advantage of zeolites over other types of adsorbents/membranes is that they typically offer good endurance to high pressure and temperature and can often tolerate harsh chemical environments. The Si-Al ratio of the zeolite framework is an important factor of the application. Zeolites with lower Si-Al ratios are more hydrophilic, while zeolite membranes with higher Si-Al ratios have fewer structural defects. High-silica zeolite membranes are greatly preferred in gas separation applications, and among them silicalite (MFI) is one of the most commonly studied.

Purification of methane from carbon dioxide (carbon dioxide is one of the major contaminants of natural gas) turns out to be an especially attractive gas separation process. Another important practical system that is involved in the process of purification of synthetic gas obtained from steam reforming of natural gas is the methane-hydrogen system. The available literature data, however, are scarce for all the systems that combine methane-hydrogen gas mixtures with the use of zeolite membranes [3, 4].

The behavior of material systems, characterized at the atomic level, can be effectively studied using classical molecular simulations. The atomistic simulation results can explain or, in some cases, substitute the experimental results. Nowadays, these simulations play an important role in the description of processes that occur in crystalline adsorbents and membranes [5, 6, 7, 8, 9, 10, 11, 12, 13]. To simulate steady-state transport of molecules, we need a method that satisfies two criteria: it must ensure the real dynamics of a system, and it must maintain the steady-state driving force at the microscopic level. The dynamics of a system can be investigated microscopically by means of molecular dynamics (MD), dynamic Monte Carlo [8] and other direct or indirect simulation methods (e.g., a very fast NP+LEMC method [11]). To preserve a constant driving force, these methods are often linked with the other techniques, e.g., dual control volume or local equilibrium Monte Carlo (LEMC) techniques [12]. There are many composite methods having many limits of applicability and numerous advantages and disadvantages, such as the gradient relaxation MD [13] or other earlier gradient techniques (e.g., [14, 15]), a self-adjusting plates technique [16], the external field MD [17] or its boundary-driven version [18], and the dual control volume grand canonical molecular dynamics (DCV-GCMD) [5]. In general, the main problem with simulating steady-state membrane transport at the atomic level is that we cannot simulate (at least virtually) bulk fluid phases of realistic size. Therefore, particle depletion/accumulation in the simulated bulk regions frequently occurs, if particle reinjection/removal steps are not applied. A sudden appearance and annihilation of molecules, however, can disturb the steady-state flux of transporting particles. Recently, we developed a simple atomistic simulation method for membrane transport that can maintain a driving force without significantly disturbing the previously-developed steady-state flux [19], while it properly mimics the common experimental situation in gas permeations measurements, where pressure is the main control parameter.

In this work, direct simulation results for steady-state gas transport through some relevant pure silica zeolite membranes are reported using our novel hybrid MD simulation method. In what follows, we briefly outline the applied transport simulation technique, and then specifications of the performed simulations and results for the separation of the technologically important CH4-H2 mixture (that might also be relevant in the development of engine fuels with high hydrogen content) will be presented in more detail.

2 Transport simulation method

The idea behind our method is based on the fact that in experimental gas permeation arrangements, pressure is the property that can be controlled relatively simply: the partial pressure for each component on the input side of the membrane and the total pressure on the permeate side. The essential point of our atomistic simulation scheme (called pressure-tuned, boundary driven molecular dynamics technique, PBD-MD) is that pressure is controlled by adjusting the particle numbers on the two sides of the membrane indirectly [19].

We started from the DCV-GCMD approach [5] containing two equilibrium control cells distinguished by unequal intensive thermodynamic parameters (such as chemical potential or pressure) on the two sides of the membrane, where random particle insertion/deletion moves are applied to maintain bulk phase properties in the control cells and thus to uphold the desired (constant) driving force through the membrane. However, to avoid such non-physical particle moves in the vicinity of the membrane affecting any eventual sorption layers on the membrane surfaces and interfering with earlier stabilized flux of the permeating particles, adjusting the particle numbers is restricted to zones far from the membrane region. This means that while all properties calculated and monitored in the control cells predominantly stem from the movements and interactions of particles in these regions, the artificial perturbation of the system is only present in the boundary regions of the simulation cell [18]. Our earlier test clearly showed that the particles inserted into the system have no “memories” of their initial velocities before they reach the interaction range of the membrane [19].

Here, we consider the system as being at a constant temperature $T$ with a fixed number of particles $N$ and a box volume $V$. To attain the target pressure $p$ (or partial pressures in the case of mixtures), regular perturbations in the number of particles are applied close to the boundaries of the total simulation box. In this way, as the chemical potentials of the control cells used in the traditional DCV-GCMD technique are connected to the corresponding partial pressures, the pressure can be controlled effectively in the control cells. In our trial-and-error type pressure-tuning approach, the particle insertion and deletion steps are performed randomly (i.e., randomly chosen particles and positions) at both ends of the simulation box, in the direction of the transport, far enough from the membrane region (in this respect, this technique is similar to the Particle Counting (PACO) method of Berti et al. [20]). One particle insertion ($+1$ case) or deletion ($-1$ case) is executed periodically if the following general criterion is satisfied:

[TABLE]

where $N_{\mathrm{control}\,\mathrm{cell}}$ is the number of particles in the control cell, and $p_{\mathrm{control}\,\mathrm{cell}}$ and $p_{\mathrm{target}}$ are the calculated pressure and the designated pressure, respectively, in the same control cell. Here, $p$ denotes partial pressure on the feed side or total pressure on the permeate side, and $N_{\mathrm{control}\,\mathrm{cell}}$ denotes the number of particles of the individual components on the feed side or the total number of particles on the permeate side. The attempts are always accepted, except for those random insertion steps which result in particle overlap (in such cases, the insertion attempt is repeated).

3 Simulation details

In the transport simulations, standard MD was used employing the leap-frog algorithm as an integrator with the time step of 2 fs. The Berendsen thermostat with the thermal coupling parameter of 200 fs (very weak coupling) was applied for controlling the temperature. The temperature was defined by subtracting off the streaming velocity from the $x$ component of the particle velocities, where $+x$ direction is the direction of the transport. The simulation box was confined by impenetrable, soft repulsive walls in the direction of the transport, while periodic boundary conditions were applied in the other two directions. (Note that the use of periodic boundary conditions in the direction of the transport is incompatible with the applied particle number adjustment procedure.) The geometry of this simulation box is depicted in figure 1. The boundary regions were considerably wider than the range of the repulsive walls on the two sides of the box. Random initial velocities were assigned to the inserted particles according to the Maxwell-Boltzmann distribution at the prescribed temperature. One particle insertion or deletion step was performed periodically in both boundary (wall) regions, outside the range of the repulsive walls. The length of the simulation period between two particle number modifications was taken as 2000 consecutive MD steps and $N_{\mathrm{control}\,\mathrm{cell}}$ and $p_{\mathrm{control}\,\mathrm{cell}}$ were collected as averages for these periods. The average pressure values were calculated by the virial expression. The membrane transport processes were simulated for at least 30 million MD steps, but in several cases twice as long runs were needed to collect a sufficient amount of transfer events (at least several hundred).

We also executed equilibrium adsorption simulations in the grand canonical ensemble (fixed chemical potential, volume and temperature) using the standard grand canonical Monte Carlo (GCMC) technique. The pressure/partial pressure of the adsorbate molecules in the gas phase was given indirectly by specifying the chemical potential of the component. The simulation box size was equal to the size of the investigated zeolite crystal and periodic boundary conditions were applied in all three spatial directions. The external surface adsorption, therefore, was not taken into account. The adsorption simulations were conducted for 200 million Monte Carlo steps.

We present results for the adsorption (equilibrium state) and transport (steady-state) of gases on all-silica MFI, LTA and DDR zeolites at $T=298.15$ and 398.15 K, and at $p=250$ and 500 (or 1000) kPa. These are built up from SiO4 tetrahedrons, with linear and/or zig-zag (sinusoidal) channels [21]. In the MFI zeolite [figure 2 (a), (b)], straight channels of elliptical cross-section extending along one direction are cross-linked by zig-zag channels of nearly circular cross-section extending in the other perpendicular direction (in both cases the channel diameters are about 0.7 nm). The crystal lattice parameters of its orthorhombic unit cell are $a=2.0090$ nm, $b=1.9738$ nm and $c=1.3142$ nm, with all the lattice angles being 90°. The pure silicon type of the LTA zeolite [figure 2 (c), (d)] has a cubic structure (all the lattice angles are 90°) with nearly spherical sodalite cages and straight channels with a diameter of about 1.1 nm between them. The crystal lattice parameter of its unit cell is $a=1.1919$ nm. The pure silica DDR zeolite [figure 2 (e), (f)] has a two-dimensional pore network with channels of slightly elliptical cross-section and with the characteristic channel diameter of about 0.8 nm. The crystal lattice parameters of its hexagonal unit cell are $a=1.3795$ nm, $b=1.3795$ nm, $c=4.0750$ nm, $\alpha=\beta=90$° and $\gamma=120$°.

The membrane structures were constructed according to the available crystallographic information from the IZA database [21] and taken to be defect-free. For the adsorption simulations, we built the model structures in the $1\times 1\times 2$, $2\times 2\times 2$ and $2\sin\gamma\times 2\times 1$ arrangements of the unit cells with MFI, LTA and DDR adsorbents, respectively. In the transport simulations, 1 or 2 unit-cell thick and 10–40 unit-cell wide (tall) membranes were used with the intention that the membrane sizes of the different zeolites should be as similar as possible. Two crystal orientations were realized in the membranes to place the (1) straight or (2) zig-zag channels of the zeolites in parallel with the $x$ direction [in LTA membranes the straight channels were rotated by 45° for case (2)]. Thicknesses and surface sizes of these membranes containing $5{-}8000$ Si and O atoms are summarized in table 1. The particles of the zeolite lattices were fixed at the crystallographic positions, so the frameworks were kept rigid [22]. In the steady-state transport simulations, the box size in the direction of the transport was set to 80 nm.

To describe the molecular interactions, the shifted and cut Lennard-Jones (LJ) pair potential was employed together with the Lorentz-Berthelot combining rule. The effective interaction potential for a pair of particles ($\alpha$ and $\beta$) was calculated as

[TABLE]

where

[TABLE]

is the standard 12-6 LJ potential, $u_{\mathrm{c}}^{\alpha\beta}$ is the value of the potential at the cut-off distance $r_{\mathrm{c}}^{\alpha\beta}$, and $\sigma$ and $\epsilon$ are the size and energy parameters of the LJ potential, respectively. The spherical cut-off distance was set to $r_{\mathrm{c}}^{\alpha\beta}=3.5\sigma^{\alpha\beta}$ leading to $u_{\mathrm{c}}^{\alpha\beta}=-0.00217478\epsilon^{\alpha\beta}$. All particles interacted with the walls according to the Weeks-Chandler-Anderson potential. This is also a shifted and cut LJ potential but it has a short cut-off $2^{1/6}\sigma^{\mathrm{wall}}$, and thus its interaction can only be repulsive. This particle-wall potential was chosen because it is convenient to be used in MD simulations. Its choice does not influence the results in the membrane region.

Single-site (i.e., not atomically detailed) models were used for the gas molecules, and the interaction parameters were taken from the literature [23, 24]. In addition, the same atomic interaction parameters were applied for different zeolites (see table 2): literature $\sigma$ and $\epsilon$ parameters [23, 25] of the zeolitic Si and O atoms were slightly adjusted to better reproduce the available equilibrium adsorption data for MFI.

4 Results

In the adsorption simulations, the equilibrium selectivity was defined as follows:

[TABLE]

Here, $q$ denotes the loading of the zeolite (adsorption inside the adsorbent in $\text{mol}\cdot\text{kg}^{-1}$, relative to the mass of the adsorbent). In the transport simulations, the permeation selectivity (or dynamical selectivity) was calculated from

[TABLE]

where $j$ is the component flux. For comparison, the permeation ratio ($R_{\mathrm{P}}$) was also evaluated; $R_{\mathrm{P}}$ is formally similar to $S_{\mathrm{P}}$ but here the pure component fluxes are used, and thus it can be regarded as the idealized limiting case.

First, the applied size and energy parameters of the potential models were tested with equilibrium adsorption calculations for CH4 and H2. Single-component adsorption isotherms were determined by GCMC simulations to select/verify the potential parameters for the pure gas components and the zeolite atoms. Here, we present the results obtained with the final parameters shown in table 2.

In the case of MFI, as expected, the calculated loadings with CH4 agree quite well with the experimental adsorption isotherm data [26] in the investigated temperature and pressure ranges [figure 3 (a)]. At lower pressures, a relatively good agreement with the available experimental data [27] can be also observed with H2 [figure 3 (b)]. There are some inaccessible cavities in this zeolite for the gas molecules, which turned out to be effectively inaccessible in the equilibrium simulations. For the pure silica LTA zeolite, we found only one experimental adsorption isotherm at room temperature or above, and only with CH4, and its data [28] were systematically overestimated in the simulations [figure 3 (c), (d)]. However, when we artificially prevented the creations of CH4 molecules inside the sodalite cages, an acceptable matching between simulation and measurement could be detected. In such a way we took into account the physical diffusion pathways in the zeolite, since the CH4 molecules are incapable of passing through the small windows of the sodalite cage. For the pure silica DDR zeolite, experimental data with CH4 were available at two temperatures [29] [figure 3 (e), (f)] and we did not find such data with pure H2. As it is seen, the simulated adsorption isotherms with CH4 deviate considerably from their experimental counterparts. In our grand canonical simulations, however, this zeolite contained a large quantity of occupiable pores or cavities, which have, according to the literature [30, 31], too small windows to allow molecule access, and so these simulations also overestimate the proper results (in these simulations, due to the complicated geometry, we could not preselect the regions where artificial molecule creations are forbidden). Accessible volume data from the literature make possible a rough estimation of an average (temperature, pressure and composition-independent) correction factor for the calculated adsorption loading of DDR: this amounts to 0.6 or 0.7 with CH4 or H2, respectively. Notwithstanding all these, we must acknowledge that we had to make a compromise at this point to keep the transferability of the pair potential parameters and highlight the differences in the simulated results that originate from the structural alterations between the zeolites.

Mixture adsorption was studied using GCMC simulations with an equimolar CH4-H2 gas system at temperatures of 298.15 and 398.15 K, and at the total pressures of 250 and 1000 kPa (more accurately, the CH4-H2 mixture was a nearly equimolar mixture, as we set the partial pressures of the components to be equal). The obtained equilibrium loading and selectivity values are shown in figure 4. The $S_{\mathrm{E}}$ values are between $\sim$4 and $\sim$80, so the degree of adsorption is noticeably higher for CH4 than for H2. We found the absolute loadings with H2 to be rather low (not shown), and this suggests some exclusion effect against the smaller gas molecules. The relatively larger $S_{\mathrm{E}}$ values for the model DDR zeolite are due to its higher loadings with CH4 (even if we consider the necessary corrections with the ratio of its inaccessible cavities). Generally, the calculated selectivity decreases with an increasing temperature and pressure.

Our simulations using the PBD-MD technique delivered results for binary gas transport through the zeolite membranes with the CH4 and H2 gases at temperatures of 298.15 and 398.15 K, and at the total feed side pressures of 250 and 1000 kPa. The permeate side pressure was always equal to 100 kPa (approximately ambient condition). As a reference, single-component transport simulations were also performed, where the feed side pressure was equal to either the partial pressure of the component or to the total feed side pressure of the corresponding mixture simulation.

Before we proceed with the permeation results, we give a validation of our PBD-MD approach for the present conditions. The concentration profile of the transporting particles is a good indicator whether the control cells can be considered as the ones representing bulk phases or not. For a few simulations with the MFI zeolite, figure 5 illustrates that there are practically constant concentrations along the $x$ direction in the control cells. Slight deviations from the horizontal line (and thus bulk character) can only be detected near the boundary region of the permeate side. In these diluted gases, the use of identical pressures or partial pressures on the feed side gives rise to very similar concentration profiles of the unlike components. The obtained concentration ratio of systems of different temperatures also properly reflects the applied temperature ratio [figure 5 (c), (d)]. It should be noted, furthermore, that the calculated total pressure on the permeate side was always correct in the simulations (the target pressures were generally reproduced within 1%), even if the partial pressures and component concentrations looked different here [e.g., figure 5 (c)]. This corresponds to the selectivity-related properties of the membrane.

The calculated permeation selectivities and permeation ratios are depicted in figure 6. It is seen that the $S_{\mathrm{P}}$ values are mostly less than 1, which means that, in their majority, the model membranes are selective for the transport of H2. Sufficiently high or low selectivity data for practical use are absent in the figure. Separation of CH4 from H2 on the permeate side is only reliable with the zig-zag model MFI membrane at a lower temperature and higher pressure, where the dynamical selectivity is almost equal to 2. On the other hand, separation of H2 from CH4 is preferred with the model DDR membranes, especially at a higher temperature and pressure ($S_{\mathrm{P}}<0.25$). In this respect, the model LTA membranes are of intermediate type. The calculated $S_{\mathrm{P}}$ values are consistently larger than the corresponding $R_{\mathrm{P}}$ values, indicating the presence of momentum coupling between the unlike gas molecules during the mixture transport and, more specifically, the effect of the prevalent adsorption of CH4. The influence on each other’s transport flux can also be recognized from the absolute values of the component fluxes because, generally, we obtained somewhat larger values for the pure component transport. Similarly to the mixture transport data (see table 3), the pure component fluxes (not shown) are the lowest for the model DDR membranes (the deviation is more significant for CH4). For the investigated membranes and conditions, the most remarkable fact is that the permeation selectivity values are much lower than the equilibrium selectivity values. It means that the adsorption preference for CH4 observed in all equilibrium cases does not predictably go hand in hand with transport selectivity in favor of CH4. This remains true even if we know that the feed side pressure in the dynamic simulations cannot be unequivocally compared to the equilibrium pressure of the adsorption simulations. In the light of the pressure gradient applied in the transport simulations across the membranes, however, it is somewhat logical that the equilibrium adsorption loadings are generally larger than the dynamic loadings of the corresponding membranes with $p_{\mathrm{feed}}=p_{\mathrm{adsorption}}$ (compare data of table 3).

The $S_{\mathrm{P}}$ values predominantly exhibit some decrease with an increasing temperature, which is attributed to a smaller decrease of the flux of H2 with temperature (cf. table 3). At the same time, the general declining pressure dependence of $S_{\mathrm{P}}$ observed here is due to a more intense increase of the flux of H2 with pressure. The only exception to this latter tendency is the transport through the model MFI membranes (with both straight and zig-zag channels) at a lower temperature, where the growth in the flux of CH4 is larger. This behavior cannot be explained by anomalous loadings of the MFI membrane in dynamic situations. As data of table 3 show, the CH4 to H2 loading ratios of the membranes in the dynamic simulations roughly follow the trends of the $S_{\mathrm{E}}$ values obtained from the equilibrium adsorption simulations. However, the much milder temperature and pressure dependences of the $R_{\mathrm{P}}$ values than those of the $S_{\mathrm{P}}$ values (see figure 6) underline the importance of the momentum coupling between different components connected to some adsorption preferences.

5 Discussion and conclusions

We used the novel PBD-MD simulation technique to study the stationary membrane transport at the molecular level. The method makes possible an accurate determination of steady-state fluxes of gases through microporous membranes. We investigated the permeation of CH4 and H2 gases across pure silica zeolite membranes (MFI, LTA and DDR types) at 25°C and 125°C and at 2.5 and 5 (or 10) bar.

The applied shifted and cut LJ potential is a less detailed interaction model for the simulations, and certainly, we could have better reproduced the available experimental data for the adsorption of pure CH4 and H2 gases by using more realistic all-atom potential models involving sophisticated inter- and intramolecular interactions. By this simple modelling, we rather preserved the transferability of the zeolite parameters and stressed the structural differences of the investigated zeolites. On the other hand, in such a way the total computational need of the executed transport simulations could be kept at a relatively low level (a couple of months with some dozens of CPU cores in a cluster computer).

We did not find experimental data for the investigated transport systems. At the same time, there would be many difficulties with the comparison of experimental and such atomistic simulation results. Experimental (synthesized) zeolite membranes typically consist of microcrystals with random orientations, and consequently, they exhibit an anisotropic pore geometry. The permeation properties of these membranes can be very different from those with oriented microcrystals [32] (these latter can be prepared in exceptional synthesis conditions). Moreover, real transport processes in zeolite membranes can occur in the intercrystalline pores. Other common problems include the presence of defects in the microcrystals and contaminations on the membrane surfaces, which may significantly affect the flux and selectivity of the membranes in experiments. Furthermore, an experimental zeolite membrane is built up from a thin zeolite layer and a macroporous supporting layer (this provides the membrane with the requisite mechanical strength); the support resistance may also influence the separation performance. From the simulation viewpoint, the main obstacle is that the extremely high demand of computational power today prevents the atomistic simulation of systems with real-life size (or nearly real-life size) membrane thicknesses (which is in the order of a micrometer). According to our former experience, the dependence of the transport properties on membrane thickness calculated with commonly applied simulation system sizes (with $1,2,3,\ldots,$ unit-cell thick membranes) cannot be easily used to forecast the magnitude of the scale-up effect, since this dependence is not necessarily linear at the atomic scale. For these reasons, we consider our present calculations only as one step towards model studies of the interplay of adsorption and diffusion in the permeation processes with realistic microporous membrane structures.

The results for the investigated systems and conditions have two striking features. First, there are no profound differences between the three microporous zeolite structures of similar channel diameter concerning both their equilibrium and transport selectivities. Second, the most important trend that is revealed from the present molecular simulations is that these model zeolites suffer very low transport selectivities for CH4 over H2 as compared to the equilibrium ones. Higher equilibrium adsorption of CH4 observed with these zeolites is more frequently linked with the permeation selectivity for H2 over CH4. The contrast between the equilibrium and transport selectivity data can be ascribed to the quite different adsorption preferences of the gas molecules, the strongly dissimilar mobilities of the components in the membrane and their dependence on the loading of the membrane. Here, a lower mobility of one gas component is compensated by a more intense adsorption capability, and the reverse situation occurs with the other component, leading to very similar calculated component fluxes, and consequently, to permeation selectivity values not far from 1. The manifestation of this compensating effect is most remarkable in the case of the model DDR zeolite, where one of the largest equilibrium loadings of the adsorbent is combined with the significantly lower fluxes of the membranes. A certain role of the different degree of temporary sorption of molecules in dynamic conditions is evidenced by the fact that for each zeolite, somewhat higher permeation selectivities were obtained for the zig-zag channels of longer transfer route than for the straight channels, which have an otherwise similar diameter. Note that the applied pair of the LTA model membranes is a special case, where the transfer route of one of these membranes is artificially lengthened through rotating the straight channels by 45°. At the same time, the latter observation with the LTA membranes suggests that the use of larger membrane thicknesses would probably result in the increase of permeation selectivity. We obtained a truly higher permeation rate of CH4 in only one case: for the model MFI membrane with zig-zag channels in the direction of the transport at a lower temperature and higher pressure. This prevalent transport of CH4 through the membrane seems to be a favorable interplay between a relatively high equilibrium selectivity and a smaller channel diameter with an optimal channel structure/orientation, which latter impedes the diffusion of H2 molecules to a greater extent. In the present systems, the transport through the membranes was more often easier for the far faster moving and smaller H2 molecules. To this end, however, there was one condition to be fulfilled: the absence of a significant steric hindrance from the adsorbed or slowly moving CH4 molecules. In our case, apparently, the mobility of the H2 molecules was not noticeably reduced due to such a blocking effect.

Acknowledgements

Present article was published within the framework of the project GINOP-2.3.2-15-2016-00053 [“Development of engine fuels with high hydrogen content in their molecular structures (contribution to sustainable mobility)”]. We gratefully acknowledge the support of the Hungarian National Research Fund (OTKA NN113527) in the framework of ERA Chemistry, too.