Selective O₂/N₂ Separation Using Grazyne Membranes: A Computational Approach Combining Density Functional Theory and Molecular Dynamics

Abstract

The separation of oxygen (O₂) and nitrogen (N₂) from air is a process of utmost importance nowadays, as both species are vital for numerous fundamental processes essential for our development. Membranes designed for their selective molecule separation have become the materials of choice for researchers, primarily due to their ease of use. The present study proposes grazynes, 2D carbon-based materials consisting of sp and sp² C atoms, as suitable membranes for separating O₂ and N₂ from air. By combining static density functional theory (DFT) calculations with molecular dynamics (MD) simulations, we address this issue through a comprehensive examination of the thermodynamic, kinetic, and dynamic aspects of the molecular diffusions across the nano-engineered pores of grazynes. The studied grazyne structures have demonstrated the ability to physisorb both O₂ and N₂, preventing material saturation, with diffusion rates exceeding 1 s⁻¹ across a temperature range of 100–500 K. Moreover, they exhibit a selectivity of ca. 2 towards O₂ at 300 K. Indeed, MD simulations with equimolar mixtures of O₂:N₂ indicated a selectivity towards O₂ in both grazynes with ca. twice as many O₂ filtered molecules in the [1],[2]{2}-grazyne and with O₂ representing ca. 88% of the filtered gas in the [1],[2]{(0,0),2}-grazyne. [1],[2]{2}-grazyne shows higher permeability for both molecules compared to the other grazyne, with O₂ demonstrating particularly enhanced diffusion capacity across both membranes. Further MD simulations incorporating CO₂ and Ar confirm O₂ enrichment, particularly with [1],[2]{(0,0),2}-grazyne, which increased the presence of O₂ in the filtered mixture by 26% with no evidence of CO₂ molecules.

1. Introduction

The Earth’s atmosphere is made up of a number of different gases, which form a perfect mixture so that life on our planet has been able to develop for millions of years. Among all the present gases, there are two that stand out above the rest: oxygen (O₂) and nitrogen (N₂), both essential for the existence of life on Earth. N₂ represents 78% of the composition of the atmosphere [1], and, without it, life as we know it could not be understood due to the so-called nitrogen cycle, where N₂ moves between the atmosphere, soil, living beings, and water, thus nourishing plants and regulating the food chain. On the other hand, O₂ represents 21% of the gases present in the atmosphere [1], and, like N₂, is one of the most important chemical elements since it is present, directly or indirectly, in almost all known biological functions. The combination of these two compounds forms a perfect mixture for the survival of living beings but, sometimes, one needs to separate them to be able to use them in different processes of technological importance. Thus, N₂ is used, for instance, in food and beverage packaging as a colorless, inert, and odorless compound in the so-called modified atmosphere packaging (MAP) [2], which reduces the oxidation of the product as well as the growth of aerobic microbes [3]. Due to its inert nature, N₂ is commonly used in chemistry to create an inert atmosphere, especially for those reactions involving highly reactive compounds [4]. Furthermore, N₂ also plays a crucial role in processes like the Haber–Bosch (HB) process, where atmospheric N₂ is used as a reactant under high pressure and temperature to produce ammonia (NH₃), which is subsequently used, along with other compounds such as carbon dioxide (CO₂) or sulfuric acid (H₂SO₄), in the production of nitrogen fertilizers. As for O₂, it is also used in MAP [2], but its applications go far beyond this. For instance, the steel industry is the largest consumer of pure O₂, using it in the process of blowing high-carbon steel [5]. In the chemical industry, commercial O₂ or oxygen-enriched air is used in the synthesis of controlled oxidation products [5]. Additionally, O₂ has medical uses, such as in oxygen therapy, which is used to treat any condition that affects the body ability to obtain and use oxygen [6].

Even though the separation of these gases—by the so-called chemical resolution—might seem simple, this statement could not be further from truth. The high stability and concomitant low reactivity of both molecules makes it difficult to separate them. In recent years, the use of different types of membranes has become a suited solution to this problem. Different types of membranes exist, predominantly categorized into organic and inorganic families. Organic membranes are principally constituted of polymers [7], whereas inorganic membranes comprise materials such as carbon and ceramics, which are mainly made of Al, Si, Ti, or Zn [8]. This latter group of membranes has been shown to be thermally stable and chemically resistant, but with very high production costs in full-scale processes [9]. Over the recent years, carbon-based membranes have gained significant attention, particularly after the discovery of graphene [10]. Membranes composed of porous materials have emerged as candidates for molecular separation, taking advantage of the difference in diffusion energies to selectively facilitate some given species permeation. Noteworthy examples include graphene oxide (GO) membranes, which have demonstrated high selectivity in separating various gas mixtures [11]. Encouraging results have also been observed on graphene nanostructures and reduced GO (rGO) applied in biogas upgrading processes [12,13]. Additionally, nanoporous graphene has exhibited notable efficiency in the separation of gas mixtures such as H₂/N₂, owing to their exceptional pore tunability [14].

Another class of carbon-based materials exhibiting favorable results are graphynes [15,16], single-sheet C-based materials like graphene, but characterized by combining C atoms with sp- and sp²-hybridizations in an ordered fashion. These structures are proposed as gas separation membranes due to their structural versatility, enabling the nanoengineering of their pores at will. The study of graphynes has given rise to a novel structural variant named grazynes [17], which comprise graphene stripes with sp² hybridization linked by acetylenic linkages. From an electronic perspective, grazyne materials exhibit novel electronic transport that occurs perpendicular to the acetylene bonds [18]. Besides, the grazyne family is also highly tunable, similarly to graphynes, either through elongation of acetylene units or graphene stripes, or by creating vacancies on the acetylene linkers. In fact, the first investigations into the application of grazynes as separation materials have yielded promising results in biogas upgrading [19] and in the separation of N₂/CO₂ from air mixtures [20].

The present study aims at exploring, by computational methods, whether the separation of N₂/O₂ gases mixtures is effectively suitable when using adapted grazyne structures as separation membranes. This study, carried out in a comprehensive, multiscale way, encompasses the examination of thermodynamic, kinetic, and dynamic aspects demonstrating the ability of both O₂ and N₂ surfaces to physisorb, with diffusion rates exceeding 1 s⁻¹ across a temperature range of 100–500 K. Moreover, they exhibit a selectivity of ca. 2 towards O₂ at 300 K. Molecular dynamics (MD) simulations with equimolar mixtures of O₂:N₂ indicated a selectivity towards O₂ in both grazynes with ca. twice as many O₂ filtered molecules. Further MD simulations incorporating CO₂ and Ar confirm O₂ enrichment, particularly with [1],[2]{(0,0),2}-grazyne, which increased the presence of O₂ in the filtered mixture by 26%, with no evidence of CO₂ molecules.

2. Computational Details

This study is based on computational simulations employing periodic density functional theory (DFT) calculations with the Vienna ab initio simulation package (VASP) [21]. The core electron density was described by the projector augmented wave (PAW) method [22], while the valence electron density expanded on a planewave basis set with a kinetic energy cutoff of 415 eV. To address exchange-correlation (xc) effects, the Perdew–Burke–Ernzerhof (PBE) xc functional was utilized [23]. Additionally, Grimme’s D3 description of dispersive forces (PBE-D3) [24] has been added to precisely capture dispersion interactions. Isolated O₂ and N₂ molecules were fully optimized at Γ k-point within a large cubic unit cell with dimensions of 10 × 10 × 10 ų, setting electronic and ionic convergence thresholds at 10⁻⁶ eV and 10⁻⁵ eV, respectively.

The investigated grazynes, named as [1],[2]{2}-grazyne and [1],[2]{(0,0),2}-grazyne, were chosen based on their favorable outcomes in biogas upgrading [19]. These grazynes are made with one-dimensional graphene stripes linked to each other by pairs of acetylenic bonds. Both structures feature a pore characterized by two consecutive acetylenic vacancies; see Figure 1 and Figure 2. The main difference lies in the arrangement of acetylenic bonds within the pore region; the [1],[2]{2}-grazyne contains a single acetylenic row between consecutive pore units, while the [1],[2]{(0,0),2}-grazyne exhibits two adjacent acetylenic rows followed by two vacancies, as indicated by the {(0,0),2} notation. For any clarification on the nomenclature of this type of materials, we refer to the literature where it is explained in detail [17].

These grazynes were fully optimized after introducing a perpendicular vacuum region of 10 Å above and below the membrane. This configuration aimed to prevent interactions between periodically repeated layers of grazyne. A Monkhorst–Pack k-point grid of 10 × 10 × 1 was used to ensure energy convergence below the chemical accuracy threshold of 1 kcal·mol⁻¹, ca. 0.04 eV. The optimized cell parameters for each grazyne model can be found in

Table 1. Optimized cell parameters, $a$ and $b$, for [1],[2]{2}- and [1],[2]{(0,0),2}-grazyne membranes.

Grazyne	a/Å	b/Å
[1],[2]{2}	7.659	18.573
[1],[2]{(0,0),2}	10.212	18.573

. The pore sizes in both grazyne structures were estimated using different methods resulting in pore areas of 59.27 Ų using C atoms as single points, a value that reduces into 44.27 Ų when subtracting the C atomic covalent radius and 52.57 Ų considering Bader’s procedure [20].

The determination of diffusion transition states (TS) corresponding to the passage of the molecules across the pore was systematically conducted. Initially, the center of mass of both N₂ and O₂ molecules was aligned with the geometric center of the pore and positioned ca. 4 Å above the layer, and oriented either perpendicular or parallel to the grazyne plane. Throughout the optimization process, the membrane furthest carbon atom from the diffusion pore was kept frozen, while all other atoms were allowed to freely relax, to avoid membrane sheet drifting. Then, the molecule was gradually brought closer to the pore until the TS and the adsorbed state were reached. Minima and TS were characterized through frequency analysis, involving the construction and diagonalization of the Hessian matrix obtained by finite atomic displacements of 0.03 Å length.

The adsorption energy for O₂ and N₂ molecules is computed through the following expression:

$$E_{\mathrm{a}\mathrm{d}\mathrm{s}}^i=E_{\mathrm{S}/i}-{(E}_{\mathrm{S}}+E_i)$$

(1)

where $E_{\mathrm{S}/i}$ represents the energy with the i^th species adsorbed on the grazyne surface, $E_{\mathrm{S}}$ is the energy of the pristine grazyne layer, and $E_i$ is the energy of the i^th species in vacuum and in their own ground state. Furthermore, the energy barrier for the diffusion across the membrane, $E_{\mathrm{b}}^i$, is obtained through

$$E_{\mathrm{b}}^i=E_{\mathrm{T}\mathrm{S},i}-E_{\mathrm{S}/i}$$

(2)

where $E_{\mathrm{T}\mathrm{S},i}$ is the energy of the TS for the i^th molecule diffusion across the grazyne pore. Both terms, $E_{\mathrm{a}\mathrm{d}\mathrm{s}}^i$ and $E_{\mathrm{b}}^i$, were gained by adding the zero-point energy (ZPE) term, defined as

$$E_{\mathrm{Z}\mathrm{P}\mathrm{E}}={\displaystyle \frac{1}{2}}\sum _j^{\mathrm{N}\mathrm{M}\mathrm{V}}h{\upsilon }_j$$

(3)

where $h$ is the Planck’s constant, and ${\upsilon }_j$ are the normal modes of vibration (NMV). The summation runs over all the NMV in the case of minima, whereas in the case of TS structures, the imaginary frequency is not included.

In addition, the rate constants associated to the diffusion process have been computed by means of the transition state theory (TST), according to the following expression:

$$r_i={\displaystyle \frac{k_{\mathrm{B}}T}{h}}{\displaystyle \frac{q_{\mathrm{v}\mathrm{i}\mathrm{b},i}^{\ne }}{q_{\mathrm{v}\mathrm{i}\mathrm{b},i}^{\mathrm{a}\mathrm{d}\mathrm{s}}}}{e}^{-{\displaystyle \frac{E_{\mathrm{b}}^i}{k_{\mathrm{B}}T}}}$$

(4)

where $k_{\mathrm{B}}$ is the Boltzmann’s constant, T is the temperature, and $q_{\mathrm{v}\mathrm{i}\mathrm{b},i}^{\ne }$ and $q_{\mathrm{v}\mathrm{i}\mathrm{b},i}^{\mathrm{a}\mathrm{d}\mathrm{s}}$ represent the vibrational partition functions corresponding to the TS and the adsorption minimum, respectively. While TST is typically used in the context of chemical reactions, it is equally applicable to diffusion processes when interactions are made/broken, as TSs can be identified on any potential energy hypersurface, regardless of the type of process involved (molecular movement or chemical bond-breaking processes). Since both states imply a strong interaction between the molecule and the membrane, only the vibrational partition function is considered, as follows:

$$q_{\mathrm{v}\mathrm{i}\mathrm{b},i}^{\mathrm{a}\mathrm{d}\mathrm{s}}={\prod }_j^{\mathrm{N}\mathrm{M}\mathrm{V}}{\displaystyle \frac{1}{1-e^{-\left(\frac{h{v}_j}{k_{\mathrm{B}}T}\right)}}}$$

(5)

Once the rates for each species are obtained, the selectivity of O₂ relative to N₂, denoted as $S_{{\mathrm{O}}_2/{\mathrm{N}}_2}$, can be determined according to,

$$S_{{\mathrm{O}}_2/{\mathrm{N}}_2}={\displaystyle \frac{r_{{\mathrm{O}}_2}}{r_{{\mathrm{N}}_2}}}$$

(6)

Furthermore, the description of desorption rates, $r_{\mathrm{d}\mathrm{e}\mathrm{s},i}$, in grazyne structures becomes possible. To achieve this, we applied the concept of a late TS within the framework of TST as follows:

$$r_{\mathrm{d}\mathrm{e}\mathrm{s},i}=S_{\mathrm{d}\mathrm{e}\mathrm{s},i}{\nu }_{\mathrm{d}\mathrm{e}\mathrm{s},i}{e}^{{\displaystyle \frac{E_{\mathrm{a}\mathrm{d}\mathrm{s}}^i}{k_{\mathrm{B}}T}}}$$

(7)

where $S_{\mathrm{d}\mathrm{e}\mathrm{s},i}$ is the desorption coefficient, with values ranging from 0 to 1, representing the probability of an adsorbed molecule being desorbed. In this study, a desorption coefficient of one was assumed for both species. The term ${\nu }_{\mathrm{d}\mathrm{e}\mathrm{s},i}$ is defined as follows:

$${\nu }_{\mathrm{d}\mathrm{e}\mathrm{s},i}={\displaystyle \frac{k_{\mathrm{B}}T}{h}}{\displaystyle \frac{q_{\mathrm{t}\mathrm{r}\mathrm{a}\mathrm{n}\mathrm{s},i}^{\mathrm{g}\mathrm{a}\mathrm{s}}{q}_{\mathrm{r}\mathrm{o}\mathrm{t},i}^{\mathrm{g}\mathrm{a}\mathrm{s}}{q}_{\mathrm{v}\mathrm{i}\mathrm{b},i}^{\mathrm{g}\mathrm{a}\mathrm{s}}}{q_{\mathrm{v}\mathrm{i}\mathrm{b},i}^{\mathrm{a}\mathrm{d}\mathrm{s}}}}$$

(8)

where $q_{\mathrm{t}\mathrm{r}\mathrm{a}\mathrm{n}\mathrm{s},i}^{\mathrm{g}\mathrm{a}\mathrm{s}}$, $q_{\mathrm{r}\mathrm{o}\mathrm{t},i}^{\mathrm{g}\mathrm{a}\mathrm{s}}$, and $q_{\mathrm{v}\mathrm{i}\mathrm{b},i}^{\mathrm{g}\mathrm{a}\mathrm{s}}$ represent the translational, rotational, and vibrational partition functions, respectively, for the i^th species in the gas phase. The term $q_{\mathrm{v}\mathrm{i}\mathrm{b},i}^{\mathrm{a}\mathrm{d}\mathrm{s}}$ was previously introduced and refers to the vibrational partition function of the same species in the adsorbed state.

Finally, the rate for a non-activated adsorption process can be determined using the Hertz–Knudsen equation

$$r_{\mathrm{a}\mathrm{d}\mathrm{s},i}={\displaystyle \frac{S_{0,i}{p}_{i}A}{\sqrt{2\pi m_i{k}_{\mathrm{B}}T}}}$$

(9)

in which $S_{0,i}$ represents the sticking coefficient that quantifies the probability that a molecule, upon reaching the membrane, will remain adsorbed. The terms $p_i,$ $m_i$, and $A$ denote the partial pressure, molecular mass, and adsorption area associated with N₂ and O₂, respectively. We define $A$ as a·b, where these parameters are reflected in Table 1. The used sticking coefficient value is the conservative value of 0.2, consistent with values reported in similar studies [25]. The vibrational frequencies of the adsorbed N₂ and O₂ molecules and their TSs on the grazyne membranes, used in the rate calculations, are provided in Tables S1–S8 of Section S1 of the Supporting Information. These data include qualitative descriptions of the associated vibrational modes.

Additionally, classical MD simulations were carried out using the large-scale atomic/molecular massively parallel simulator (LAMMPS) package [26]. The simulation boxes used were 61.272 × 148.584 × 140 ų for the [1],[2]{2}-grazyne and 81.696 × 148.584 × 140 ų for the [1],[2]{(0,0),2}-grazyne. The grazyne membranes were located at half the z axis length along the xy plane, creating two sections with equal volumes. The adaptive intermolecular reactive empirical bond-order (AIREBO) potential was employed to define the interactions of the grazyne layer [27], thanks to its good representation of interatomic interactions within materials and compounds, predominantly constituted of carbon and hydrogen [28,29]. In the case of gas molecules, the intermolecular interactions were modelled by Lennard-Jones 12-6 potential plus a Coulomb interaction term [30], as follows:

$$V_{non-bonded}\left(r_{ij}\right)={\displaystyle \frac{q_i{q}_j}{4\pi {\epsilon }_0{r}_{ij}}}+4{\epsilon }_{ij}\left(\right({\sigma }_{ij}/r_{ij}{)}^{12}-({\sigma }_{ij}/r_{ij}{)}^6)$$

(10)

where ${\epsilon }_0$ is the vacuum permittivity, ${\epsilon }_{ij}$ and ${\sigma }_{ij}$ are the Lennard-Jones parameters of the ij interaction, $r_{ij}$ is the interatomic distance between pairs, and $q_i$ and $q_j$ denote the partial atomic charges of the ij pair species. All the values used are listed in

Table 2. Non-bonded parameters for O₂, N₂, CO₂, and Ar.

Molecule	$\mathbf{Atom},\mathit{i}$	${\mathit{q}}_{\mathit{i}}$/e	${\mathit{\epsilon }}_{\mathit{i}}·{\mathit{k}}_{\mathit{b}}^{-1}/$eV	${\mathit{\sigma }}_{\mathit{i}}$/Å
O2	O	−0.1130	49.000	3.020
	dummy	+0.2260	—	—
N2	N	−0.4820	36.400	3.320
	dummy	+0.9640	—	—
CO2	C	+0.6512	28.129	2.757
	O	−0.3256	80.507	3.033
Ar	Ar	0.0	119.800	3.405

and correspond to the transferable potentials phase equilibria (TraPPE) force field for the oxygen [31], which establishes a 1.210 Å bond length between oxygens, accompanied by a massless dummy atom positioned at the molecular center of mass. A charge of −0.113 e was assigned to each of the two oxygen atoms, balanced by the charge assigned to the dummy atom of 0.226 e. The N₂ molecules were modeled based on the methodology proposed by Murthy et al. [32], with a bond length of 1.098 Å. In order to account for the experimental quadrupole moment inherent to the molecule, a massless dummy atom was placed at the center of mass with a charge of 0.964 e, thereby equalizing the charges of −0.482 e assigned to each nitrogen. The CO₂ molecule was modeled according to the elementary physical model (EPM2) [33], in which the oxygen atoms are located at a bond length of 1.149 Å with respect to the carbon atom, resulting in a linear molecular geometry. The charge assigned to the C atom is 0.6512 e while on each oxygen atom the charge corresponds to −0.3256 e. Finally, Argon atoms have been modeled as a single particle with zero charge [30]. The unlike pair parameters were calculated using the arithmetic Lorentz–Berthelot combination rules. During the MD simulations, the molecules were treated as rigid entities.

All MD simulations were conducted under the canonical NVT ensemble, where the number of particles, N, the volume, V, and the temperature, T, of the system were fixed. Different values of pressure were studied to assess its significance in the process. To control the temperature of the system and keep it constant throughout the simulation, the Nose–Hoover thermostat was employed [34,35,36]. Initially, the gas mixture was thermalized for 100 ps, a process carried out avoiding gas–membrane interaction. To do so, the gas molecules are confined between two rigid walls, one directly above the grazyne membrane and another one at a distance ranging from 20 Å to 50 Å from the grazyne, depending on the desired pressure. Once the desired temperature was reached, the lower wall was removed and relocated at the bottom edge of the simulation box, then the diffusion process was started and was sustained over 500 ps, time that constitutes the production stage.

The number of molecules that diffuse across the surface within the production stage is commonly used to compute the membrane permeability, $P_i$, defined as [14],

$$P_i={\displaystyle \frac{N_i}{A·t}}$$

(11)

where $N_i$ are the moles of net permeated gas molecules, $A$ corresponds to the membrane area, and $t$ corresponds to the simulation time.

3. Results

3.1. DFT Results

Before presenting the results obtained from DFT calculations, a series of considerations are needed. First, to determine whether a molecule is capable of passing through the pore, it is desirable to obtain values below 1 eV for diffusion energy barriers. In this sense, previous studies have determined that values below 0.3 eV would be adequate to ensure that molecules diffuse through these materials [19]. Furthermore, it is necessary to obtain moderate $E_{\mathrm{ads}}^i$ values, as low adsorption energies would result in insufficient gas attraction, whereas high adsorption energies could lead to chemisorption processes, deviating from the desired physisorption mechanism. Chemisorption would imply that the molecules could remain captured on the grazyne membrane surface instead of permeating through, making grazyne membranes effectively scrubber materials, and to act as membranes, it would be necessary to find suitable techniques to regenerate and clean grazynes from these adsorbed molecules.

Moreover, when investigating single molecules diffusion, scenarios involving the simultaneous presence of two molecules in the same pore are excluded, given their low probability, and the repulsive interactions they would have, frustrating their diffusion through the grazyne structure. Consequently, single passage of O₂ and N₂ molecules across grazyne membranes was studied yet considering two possible molecular orientations for each molecule type, either with the molecular axis being perpendicular or parallel to the grazyne plane. The obtained energetic values are listed in

Table 3. Adsorption energies, $E_{\mathrm{a}\mathrm{d}\mathrm{s}}^i$, and diffusion energy barriers, $E_{\mathrm{b}}^i$, for O₂ and N₂ on both grazyne structures, for perpendicular and parallel orientations. All values are given in eV.

Grazyne	Perpendicular				Parallel
	$E_{\mathrm{ads}}^{O_2}$	$E_{\mathrm{b}}^{O_2}$	$E_{\mathrm{ads}}^{N_2}$	$E_{\mathrm{b}}^{N_2}$	$E_{\mathrm{ads}}^{O_2}$	$E_{\mathrm{b}}^{O_2}$	$E_{\mathrm{ads}}^{N_2}$	$E_{\mathrm{b}}^{N_2}$
[1],[2]{2}	−0.41	0.06	−0.18	0.22	−0.41	0.06	−0.17	0.27
[1],[2]{(0,0),2}	−0.32	0.02	−0.15	0.20	−0.30	0.02	−0.13	0.29

. For both structures, the adsorption of O₂ molecules is stronger than that of N₂, regardless of the considered orientation. Another positive income comes from the fact that both species physisorb ensuring no blockage of the membrane. In addition, diffusion O₂ barriers are totally non-sensitive to the molecular orientation whereas N₂ slightly prefers the perpendicular path. According to the energetic values reported in Table 3, the [1],[2]{(0,0),2}-grazyne membrane exhibits lower diffusion barriers than [1],[2]{2}-grazyne.

Using the values from Table 3, the energy profile for the diffusion of the studied molecules has been represented in Figure 3. It is evident that O₂ experiences lower energy barriers for both conformations with respect to the N₂ cases. Among the studied conformations, the perpendicular orientation appears to be more suited to carry out the diffusion through the grazyne membranes for both molecules.

The TS visual inspection helps at identifying different situations. For instance, TS (cf. Figure 4) for N₂ is obtained when the adsorbate is positioned at 0 Å from the pore center, with the layer in a completely flat conformation. In the case of O₂, the TS corresponds to the molecule slightly tilted from the completely perpendicular orientation in both structures.

The adsorbed states of the molecules also appear to be distinct. For O₂, grazyne bulges toward the molecule, indicating an attraction between the two parts. Conversely, for N₂ the same grazyne bulges away from the molecule, resulting in a repulsive interaction, as observed in Figure 5. Apparently, the energy differences shown in Table 3 could be attributed to these structural discrepancies obtained in the TS and adsorbed state of each molecule. The planar configuration observed in the TS and the bulged structure in the adsorbed state are not a result of constraints applied during the simulations. Instead, these geometries arise naturally as the system evolves to minimize its energy. This behavior reflects the response of the grazyne membranes to the interaction with gas molecules under the simulated conditions.

To further analyze the factors affecting the diffusion of N₂ and O₂ through the grazyne membranes, the electrostatic potential maps for both grazyne structures were acquired and shown in Figure 6. There, it is evident that the electrostatic potential is uniformly distributed around the atoms of the material, and nearly zero within the pore region, implying that the diffusion is not biased by the electrostatic potential, yet in turn, primarily governed by the pore size and the size of the molecule.

3.2. Rate Constants and Selectivity

According to the energy barriers obtained in the previous section, O₂ should pass more easily than N₂ on both grazyne membranes. In order to further verify this hypothesis, rate constants have been calculated for each species at different temperatures by means of Equation (4) and represented in Figure 7. For each rate constant the two orientations studied—perpendicular and parallel—have been considered. Clearly, O₂ exhibits a higher diffusion in both grazynes, reaching values ca. 12 s⁻¹ at 500 K. For both O₂ and N₂, the [1],[2]{2}-grazyne yields the largest rate constant values, even though the difference with the values obtained in [1],[2]{(0,0),2}-grazyne is very small. Furthermore, with an increase in temperature, the rate constant also rises, as expected. The temperature trend observed can be explained from a physical perspective, where higher temperatures provide molecules with greater kinetic energy and velocity, thereby facilitating the diffusion across the pores. According to the different rate constants, it is reasonable to consider both membranes for O₂/N₂ separation and identify the appropriate range of temperatures and pressures that enhance the separation.

Figure 8 represents the selectivity of both membranes towards O₂ as a function of temperature. The results show that, as the temperature increases, both membranes become less selective. Thus, at 150 K, the selectivity value is 3.5 and 3.2 for [1],[2]{(0,0),2}-grazyne and [1],[2]{2}-grazyne, respectively. However, at 300 K, both grazyne membranes achieve a selectivity value of 1.8, indicating that approximately twice as much O₂ would be filtered from the mixture, ensuring the capability of purifying a mixture of these gases even at room temperature.

3.3. Kinetic Phase Diagrams

The analysis of the ability of materials to adsorb molecules is essential in designing effective filtration systems. The effectiveness of a material in this regard depends on its ability to adsorb target molecules upon, and if so, in the regeneration for a continued use. To evaluate this process, we created Kinetic Phase Diagrams (KPDs) that illustrate the adsorption and desorption tendency of N₂ and O₂ on both grazynes as a function of the temperature and pressure of the system [37,38]. Under equilibrium conditions, the rates of adsorption and desorption are equal. This equilibrium prompts the search for specific temperatures and pressures that achieve this balance. Figure 9 clearly shows two different regions; one where molecular adsorption is preferred, and another region where desorption is preferred, separated by the equilibrium line. Keeping that in mind, it is worth noting that O₂ exhibits a superior adsorption range on both types of grazynes, consistent with the stronger $E_{\mathrm{ads}}^i$ values reported in Table 3, vide supra.

The [1],[2]{(0,0),2}-grazyne displays lower adsorption range, positioning it as a slightly more favorable option for operating at high pressure and temperature conditions. However, what is most significant is that under standard conditions (p = 1 atm and T = 300 K), neither of the molecular species exhibit a propensity for adsorption, making these grazynes potential candidates for their use as separation membranes, as they would not, in principle, get their pores poisoned by gas adsorption.

3.4. Molecular Dynamics

Molecular dynamics (MD) simulations allow a direct comparison with results obtained by DFT calculations and check whether both methodologies are in line and complement each other. While MD focuses on capturing the dynamic aspects of processes, DFT has so far been limited to addressing thermodynamics and static kinetic approximations. The results of these MD simulations—the number of molecules of each species crossing the membrane—have been used to calculate the permeability of grazynes from Equation (11). Here, two types of simulations can be differentiated depending on the gas mixture that has been used. In the first type, a gas mixture with equimolar composition of O₂ and N₂ has been modeled (i.e., 64 molecules of O₂ and 64 molecules of N₂) and the diffusion across both membranes was evaluated at different pressures. For the first set of simulations, three different simulations were performed for each pressure, with the variation among replicates arising from differences in their initial states. These additional simulations serve as independent replicates, enabling assessment of the robustness of the results and confirming that the observed trends are not merely due to probabilistic situations. The results of these simulations are shown in

Table 4. Number of permeated gas molecules in [1],[2]{2}-grazyne for an equimolar initial mixture of O₂/N₂ at T = 300 K and at different pressures, for the three different MD runs.

	[1],[2]{2}-grazyne
P (atm)	30.6			20.1			14.9			11.9
#R	1	2	3	1	2	3	1	2	3	1	2	3
#O2	37	40	34	28	30	32	32	24	29	28	14	18
#N2	17	25	21	14	14	13	11	13	13	6	9	9
%O2	68.5	61.5	61.8	66.7	68.2	71.1	74.4	64.9	69.0	82.4	60.9	66.7

and

Table 5. Number of permeated gas molecules in [1],[2]{(0,0),2}-grazyne for an equimolar initial mixture of O₂/N₂ at T = 300 K and at different pressures, for the three different MD runs.

	[1],[2]{(0,0),2}-grazyne
P (atm)	31.0			20.2			15.0			11.8
#R	1	2	3	1	2	3	1	2	3	1	2	3
#O2	25	33	20	29	26	22	21	13	13	14	15	20
#N2	5	12	12	4	9	5	5	5	4	6	2	1
%O2	83.3	73.3	62.6	87.8	74.3	81.5	80.8	72.2	76.5	70.0	88.2	95.2

. Here, O₂ exhibited a greater capacity for crossing the membrane, a fact that is consistent with a lower DFT energy barrier associated with O₂ (cf. Table 3). To be exact, approximately twice as many O₂ molecules were filtered in [1],[2]{2}-grazyne. In the [1],[2]{(0,0),2}-grazyne, the difference is even more noticeable where O₂ represented up to ca. 95% of the total filtered gas species, which implies an increase of ca. 45%. These differences are maintained for all pressure values, where a small tendency to decrease the number of filtered molecules with decreasing pressure can be seen. In fact, the values reported in Table 4 show that by decreasing the operative pressure, the membrane generates a purer O₂ phase using the [1],[2]{2}-grazyne. Note that the pressure indicated for each simulation represents only the initial pressure of the gas on the feeding side. As seen in the MD (see Trajectory S1 in the Supporting Information), the deformation and vibration of the membrane over time slightly modifies the cell volume and shape at each timestep, making it challenging to accurately evaluate the feed pressure throughout the simulation. When one focuses on the replicates, for each pressure, oscillations are observed in the number of filtered molecules. However, a clear trend is observed with O₂ consistently exhibiting the highest filtration rate.

Tracking the molecules that come across the membrane at each timestep enables the generation of Figure 10, where the number of molecules in the permeated side is represented over time. The thermalization process, occurring from 0 to 100 ps, is omitted. Initially, both molecules exhibit a comparable diffusion rate. However, before reaching 200 ps, there is a noticeable shift in trends. Nitrogen diffusion appears to decelerate, whereas oxygen maintains its initial diffusion pattern unaltered.

By using the data reported in Table 4, the permeability was computed and represented in Figure 11. Comparing both grazynes, it becomes evident that [1],[2]{2}-grazyne exhibits higher permeability values for both molecules, thus being a less restrictive membrane to O₂/N₂ diffusion. Furthermore, O₂ molecules displayed elevated permeability values across both grazynes, underscoring their superior diffusion capacity.

3.5. O₂ Direct Air Capture

Considering the results obtained through DFT and MD analyses, it appears plausible to achieve the separation of these molecules using the studied grazynes, thanks to their different diffusion rates and permeability, with oxygen showing greater ease of passage through the membranes. However, the study carried out up to this point only considers N₂ and O₂ molecules, when there are more compounds in the air that can affect the separation. To have a more representative view regarding the separation of these molecules, a new set of MD simulations were carried out considering, apart from N₂ and O₂, CO₂ molecules and Ar atoms—the main compounds after N₂ and O₂ in air—plus the mixture was made so that the composition was as close as possible to that of atmospheric dry air.

Table 6. Gas mixture composition used in the MD simulations.

Molecule Type	Number of Molecules	%
N2	7800	77.97
O2	2100	20.99
Ar	100	1.00
CO2	4	0.04

shows the number of molecules and their contribution to the four components initial mixture.

The size of the simulation box was expanded to accommodate all the molecules, being of 152.524 × 373.576 × 2996.384 ų for the [1],[2]{2}-grazyne and 204.240 × 373.575 × 2996.384 ų for the [1],[2]{(0,0),2}-grazyne.

Table 7. Number of molecules filtered through the grazyne membranes for a N₂/O₂/Ar/CO₂ mixture.

	[1],[2]{2}-grazyne		[1],[2]{(0,0),2}-grazyne
	#Permeated Molecules	%Filtered Gas	#Permeated Molecules	%Filtered Gas
N2	562	60.30	231	51.11
O2	351	37.66	213	47.12
Ar	18	1.93	8	1.77
CO2	1	0.11	0	0.00

brings together the results obtained in these simulations. The trend observed in the previous simulations is replicated also here, where the filtered gas was notably enriched in O₂. Particularly noteworthy is the result for the [1],[2]{(0,0),2}-grazyne, where oxygen increases its presence up to 47% when it only constituted 21% in the initial mixture. On the other hand, N2 presence is reduced to 60% and 51% in the [1],[2]{2}- and [1],[2]{(0,0),2}-grazyne, respectively. However, a downside of these simulations is the permeation of CO₂, which, despite starting with only four initial molecules, managed to penetrate the [1],[2]{2} membrane, as anticipated based on previous findings where CO₂ diffusion was documented [19,20]. Moreover, the quantity of filtered molecules is greatly influenced by their interaction with the membrane, particularly when the approach occurs perpendicular or parallel to its surface. This fact underscores the significance of understanding the molecular dynamics and structural characteristics of the membrane to predict and optimize separation processes effectively. The [1],[2]{(0,0),2}-grazyne remained impermeable to CO₂ initially, though with prolonged simulation, some CO₂ molecules might penetrate this membrane, akin to the other model. While argon atoms also crossed the membranes, their inert nature renders them non-threatening. From a broader perspective, [1],[2]{2}-grazyne allows a greater number of total molecules to pass compared to its counterpart, indicating its less restrictive nature. Specifically, 932 molecules were filtered in [1],[2]{2}, whereas only 452 were filtered in [1],[2]{(0,0)2}. Despite this, the gas filtered through [1],[2]{2} is less enriched in O₂, representing only 37%. According to these results, the effectiveness of both membranes to generate richer O₂ mixtures directly from N₂/O₂ and dry air mixtures is notorious.

In this work, we focused on anhydrous conditions for the gas separation study, as dry gases are commonly used in gas separation tests. Consequently, the effects of moisture or air vapors were not explicitly considered. However, we acknowledge the importance of considering these factors, as seen, e.g., on water on the adsorption and diffusion processes through graphene [39]. Future advances in the present work should thus consider the influence of moisture, vapors, and various impurities on the performance of grazyne-based membranes, providing a more comprehensive understanding of their behavior under more realistic conditions.

4. Conclusions

The present DFT and MD calculations analyzing N₂ and O₂ permeation through nano-engineered grazyne membranes in a holistic fashion, considering thermodynamic, kinetic, and dynamic aspects yield consistent results, indicating the potential use of grazynes to separate N₂ and O₂ from air. DFT calculations show the physisorption of both molecules on the grazyne membrane models, preventing the pore obstruction by chemical adsorption, further confirmed by KPD under working temperature and pressure conditions. Additionally, the DFT calculations suggest that the two types of molecules could permeate the membrane in two different orientations: perpendicular and parallel. In this sense, O₂ exhibited a significantly higher diffusion capacity compared to N₂, attributed to its smaller diffusion barriers. This trend is reflected in the diffusion rates where O₂ obtained systematically higher rates than N₂ but with both rates higher than 1 s⁻¹. Despite the significant differences in the values, working at low temperatures would be preferable to enhance selectivity towards O₂. Finally, the electrostatic potential maps indicate that the diffusion of N₂ and O₂ through the grazyne membranes is primarily influenced by the pore size and the size of the molecule rather than by electrostatic interactions between grazyne and molecules’ electron densities.

MD simulations with equimolar mixtures of O₂:N₂ clearly indicate a selectivity towards O₂ in both grazynes. To be precise, approximately twice as many O₂ molecules were filtered in [1],[2]{2}-grazyne whereas in the [1],[2]{(0,0),2}-grazyne, the difference is even more pronounced, with O₂ representing up to ca. 95% of the filtered gas, implying an increase of ca. 45%. These data align perfectly with the selectivity obtained through DFT calculations, where a selectivity towards O₂ close to 2 at 300 K was observed. The simulations also show promising results in the case of dry air. Here, the trend observed for the O₂/N₂ mixture was confirmed, with O₂ being enriched in the permeate gas for both grazynes. Particularly noteworthy is the result for the [1],[2]{(0,0),2}-grazyne, where oxygen increases its presence by 26%, nearly constituting half of the filtered gas (i.e., 47%), despite only constituting 21% of the initial gas. On the other hand, [1],[2]{2}-grazyne allows a greater number of molecules to diffuse across the membrane, indicating its less restrictive nature. However, despite this difference, the gas filtered through [1],[2]{2} is less enriched in O₂, representing only 37%. Nonetheless, the objective of obtaining enriched oxygen mixtures remains achievable.

All in all, the present study, which explores the energetics, kinetics, and dynamics of O₂ separation from N₂ using engineered grazyne membranes, paves the way for their application. However, the results underscore that it is essential to control the operating conditions to maximize the separation efficiency. Further research may be necessary to understand the separation behavior and consider the possible impact of other species present in untreated O₂/N₂ streams.