Multilayer Graphtriyne Membranes for Separation and Storage of CO2: Molecular Dynamics Simulations of Post-Combustion Model Mixtures

The ability to remove carbon dioxide from gaseous mixtures is a necessary step toward the reduction of greenhouse gas emissions. As a contribution to this field of research, we performed a molecular dynamics study assessing the separation and adsorption properties of multi-layered graphtriyne membranes on gaseous mixtures of CO2, N2, and H2O. These mixtures closely resemble post-combustion gaseous products and are, therefore, suitable prototypes with which to model possible technological applications in the field of CO2 removal methodologies. The molecular dynamics simulations rely on a fairly accurate description of involved force fields, providing reliable predictions of selectivity and adsorption coefficients. The characterization of the interplay between molecules and membrane structure also permitted us to elucidate the adsorption and crossing processes at an atomistic level of detail. The work is intended as a continuation and a strong enhancement of the modeling research and characterization of such materials as molecular sieves for CO2 storage and removal.


Introduction
The persistent growth of greenhouse gas emissions, the debate about the role of anthropic activities in connection with climate change [1], and the global warming phenomenon [2] have stimulated the search for "clean" technologies that can achieve carbon dioxide removal from gas mixtures such as the flue gases generated after combustion. Selective adsorption using porous materials is a promising way to capture CO 2 , which is mostly generated by fossil fuel combustion, in order to mitigate the greenhouse effects associated with its excessive concentration in the atmosphere [3][4][5]. This method is favored in terms of its simplicity and lower implementation costs, in comparison with the more traditional aqueous chemical absorption [6].
In practice, it is a very difficult or even impossible task to experimentally synthesize, characterize, and evaluate the performance in terms of the CO 2 capture and separation of all possible materials. Up to this point in the process, computer modeling and simulations play an important role in material design and development, prior to the experimental stage [29,30]. Due to the variety of interactions between molecules and materials, generic force fields like AMBER [31] and UFF [32] often poorly model the particular system. Thus, parts of the potential energy function must be developed or refined on purpose, using the available theoretical and experimental data. Recently, various force fields specific to graphene and its derivatives have been developed [33,34], as well as those for MOFs [35][36][37], zeolites [38], and other polymeric materials [39], to identify molecular interactions providing realistic predictions of relevant adsorption dynamics and the transport properties of gas under consideration.
However, intermolecular interaction parameterization is a delicate task and the accurate formulation of force fields is an indispensable work: they must be reliable in terms of the full space of the relative configurations of involved partners and must be expressed in a proper analytical form, to permit extensive molecular dynamics (MD) simulations under a variety of conditions of interest. The complete achievement of this knowledge is a very critical question, one that is also difficult to answer for relatively simple systems. This has been one important target of our recent research activity. We have found that the strength, radial, and angular dependencies of the basic interaction components are definable by semi-empirical and empirical functions, the parameters of which relate to the fundamental chemical-physical properties of the interacting partners. Under these conditions, such functions become scaling laws, the involved parameters assume a defined transferability degree, and both gas-gas and gas-layer intermolecular interactions become representable in an internally consistent way.
Permeability and selectivity are the two main aspects by which to determine whether a membrane can be effective for gas separation. It is already well known that permeability is inversely proportional to the thickness of the membrane. Therefore, a single-atom-thick planar membrane may have great potential for gas separation if tailored to be selective for a given molecule [40,41]. The γ-graphynes are single atomic layers belonging to the class of carbon allotropes wherein hexagonal carbon rings are connected by carbon chains containing a variable number of C-C acetylenic bonds. The γ-graphynes exhibit similar properties to graphene, but the pores are uniformly distributed and have adjustable dimensions [42]. Moreover, graphynes have lower dispersion forces that minimize aggregate formation among the layers and molecules. The synthesis and characterization techniques of graphynes have actively been developed over the last few years [43][44][45][46][47][48]. In our previous works [49,50], force fields related to gas adsorption on γ-graphynes have also been developed and tested using accurate ab initio calculations. Therefore, in this work, we will perform extended molecular dynamics (MD) simulations on a wide range of conditions that are typical of post-combustion gaseous mixtures containing CO 2 /N 2 /H 2 O, in order to characterize the separation properties of multi-layer graphtriyne membranes, a γ-graphyne characterized with three consecutive C-C acetylenic bonds.

Methods
The classical MD simulations of graphtriyne and gas mixtures were performed by enforcing periodic boundary conditions with a simulation box with the dimensions  Three different types of arrangements, with 1, 2, and 3 membrane layers at  four different temperatures (333, 353, 373, and 400 K), were the subject of simulations with CO 2 /N 2 /H 2 O gaseous mixtures, with equal concentrations of CO 2, N 2 , and H 2 O. The structures of the graphtriyne membranes were taken from Ref. [49], where they had been optimized using periodic DFT calculations (the structural details are reported in the Supplementary Information, in Figure S1).
As previously reported by the authors of [50], the intermolecular potential has been defined as a combination of electrostatic and non-electrostatic components. The first component is represented by the sum of Coulomb interactions between the partial point charges located on each molecular frame, the anisotropic distribution of which accounts for the permanent electric multipole of each partner. In particular, for H 2 O, the representation reported by the authors of [51] was adopted, with a charge distribution correctly reproducing the dipole moment of water in the gas phase (1.85 D) [52], while a three-charge-site N 2 model [53] and a five-charge-site CO 2 model [54] were used (see Figure 1).
The classical MD simulations of graphtriyne and gas mixtures were perfor enforcing periodic boundary conditions with a simulation box with the dimension Å × 62.523 Å × 280.0 Å. The box contained graphtriyne membrane(s) with the dim 72.210 Å × 62.523 Å, placed at the midway point perpendicularly to the z-direction different types of arrangements, with 1, 2, and 3 membrane layers at four d temperatures (333, 353, 373, and 400 K), were the subject of simulations with CO2/ gaseous mixtures, with equal concentrations of CO2, N2, and H2O. The structure graphtriyne membranes were taken from Ref. [49], where they had been optimize periodic DFT calculations (the structural details are reported in the Supplem Information, in Figure S1).
As previously reported by the authors of [50], the intermolecular potential h defined as a combination of electrostatic and non-electrostatic components. T component is represented by the sum of Coulomb interactions between the parti charges located on each molecular frame, the anisotropic distribution of which a for the permanent electric multipole of each partner. In particular, for H representation reported by the authors of [51] was adopted, with a charge distr correctly reproducing the dipole moment of water in the gas phase (1.85 D) [52], three-charge-site N2 model [53] and a five-charge-site CO2 model [54] were us Figure 1).

Figure 1.
A schematic view of carbon dioxide, nitrogen, and water molecules, with th charges (a.u.), angles, and bond lengths adopted. Charge-sites (Q) are represented by the gray spheres.
The second (non-electrostatic) component of the intermolecular potential between pairs of gas molecules and between the gas molecules and the membr determined by the balance of size repulsion with dispersion and induction attrac previously, it has been described by pair-wise additive contributions, arising fr different interaction centers distributed again on each molecular frame, and form using the improved Lennard-Jones (ILJ) potential function [55]. The ILJ exp adopted in place of the much simpler Lennard-Jones function (usual in gener fields) describes non-electrostatic intermolecular interactions in a far more accur [56][57][58]. All the ILJ parameters used in this work, predicted from the polari component assigned at each interaction center, were tested and fine-tuned us available experimental findings, exploiting the comparison of interaction energi the results of high-level ab initio calculations [49,50]. A cut-off distance for the electrostatic interactions was set to be equal to 15 Å. Given the periodic bo conditions, the smoothed particle mesh Ewald method, as implemented in the DL software (see below), was applied to accurately account for long-range elect interactions [59].
All MD simulations were performed using the DL_POLY molecular dy program [60] in the canonical (NVT) ensemble, employing the Nosé-Hoover the with periodic boundary conditions in all directions. Each simulation was carried o period of 5.5 ns, after a 0.5 ns equilibration period, with a fixed time step of 1 trajectory data for the statistics were collected every 2 ps. Seven different amount have been loaded into the simulation box for every system to characterize the influ The second (non-electrostatic) component of the intermolecular potential, acting between pairs of gas molecules and between the gas molecules and the membranes, is determined by the balance of size repulsion with dispersion and induction attraction. As previously, it has been described by pair-wise additive contributions, arising from the different interaction centers distributed again on each molecular frame, and formulated using the improved Lennard-Jones (ILJ) potential function [55]. The ILJ expression, adopted in place of the much simpler Lennard-Jones function (usual in generic force fields) describes nonelectrostatic intermolecular interactions in a far more accurate way [56][57][58]. All the ILJ parameters used in this work, predicted from the polarizability component assigned at each interaction center, were tested and fine-tuned using the available experimental findings, exploiting the comparison of interaction energies with the results of high-level ab initio calculations [49,50]. A cut-off distance for the ILJ and electrostatic interactions was set to be equal to 15 Å. Given the periodic boundary conditions, the smoothed particle mesh Ewald method, as implemented in the DL_POLY software (see below), was applied to accurately account for long-range electrostatic interactions [59].
All MD simulations were performed using the DL_POLY molecular dynamics program [60] in the canonical (NVT) ensemble, employing the Nosé-Hoover thermostat with periodic boundary conditions in all directions. Each simulation was carried out for a period of 5.5 ns, after a 0.5 ns equilibration period, with a fixed time step of 1 fs; the trajectory data for the statistics were collected every 2 ps. Seven different amounts of gas have been loaded into the simulation box for every system to characterize the influence of pressure on the observable objective of the present investigation. The resulting gas pressure was computed using the Peng-Robinson equation of state [61]. In order to mimic post-combustion conditions, the initial pressure values were chosen to be lower than 5.5 atm (the details are reported in the Supplementary Information, in Tables S1-S4). At the beginning of the simulations, equal amounts of the gas molecules were randomly distributed into the two regions of the box, in such a way that there was no pressure gradient inside the box ( Figure 2). The membranes were considered a frozen framework and the gas molecules were treated as rigid bodies. The gas molecules could cross the membrane multiple times in both directions of the z-axis during the simulation. The number of permeation events was then monitored, along with the z-density and radial distribution function profiles. All graphical representations of the molecular trajectories were processed using the VMD package [62].
pressure on the observable objective of the present investigation. The resulting gas pressure was computed using the Peng-Robinson equation of state [61]. In order to mimic post-combustion conditions, the initial pressure values were chosen to be lower than 5.5 atm (the details are reported in the Supplementary Information, in Tables S1-S4). At the beginning of the simulations, equal amounts of the gas molecules were randomly distributed into the two regions of the box, in such a way that there was no pressure gradient inside the box ( Figure 2). The membranes were considered a frozen framework and the gas molecules were treated as rigid bodies. The gas molecules could cross the membrane multiple times in both directions of the 007A-axis during the simulation. The number of permeation events was then monitored, along with the -density and radial distribution function profiles. All graphical representations of the molecular trajectories were processed using the VMD package [62].

Gas Permeability
The simulations were considered to be initialized after the equilibration steps when the permeation events that occurred were monitored and counted. The numbers of permeation events were then plotted against the simulation time. By definition, the slope of such plots represents an estimation of the gas permeation rate, measured in units of molecules ps −1 . Using these data, the gas permeance was then calculated by dividing the permeation rate by the corresponding pressure and by the area of the membrane. The gas permeances of the single layer are reported in Figure 3, using the gas permeance unit (GPU), where 1 GPU is equal to 3.35 × 10 −10 mol m −2 s −1 Pa −1 . Figure 3A shows that the N2 permeances are not affected by the pressures. Meanwhile, the CO2 and H2O permeances are higher at low pressure and show a relatively flat trend at increasing pressures. However, for all gas molecules, the permeances do not vary much at pressures higher than 2 atm.

Gas Permeability
The simulations were considered to be initialized after the equilibration steps when the permeation events that occurred were monitored and counted. The numbers of permeation events were then plotted against the simulation time. By definition, the slope of such plots represents an estimation of the gas permeation rate, measured in units of molecules ps −1 . Using these data, the gas permeance was then calculated by dividing the permeation rate by the corresponding pressure and by the area of the membrane. The gas permeances of the single layer are reported in Figure 3, using the gas permeance unit (GPU), where 1 GPU is equal to 3.35 × 10 −10 mol m −2 s −1 Pa −1 . Figure 3A shows that the N 2 permeances are not affected by the pressures. Meanwhile, the CO 2 and H 2 O permeances are higher at low pressure and show a relatively flat trend at increasing pressures. However, for all gas molecules, the permeances do not vary much at pressures higher than 2 atm. However, upon closer scrutiny, the data in Figure 3A, although shown in the limited pressure range considered, seem to suggest a minimum level of permeance, followed by a rebound profile. Such a trend, in principle, is to be expected, due to the increase in the frequency of collisions that is typical of increasing pressures. However, more simulations However, upon closer scrutiny, the data in Figure 3A, although shown in the limited pressure range considered, seem to suggest a minimum level of permeance, followed by a rebound profile. Such a trend, in principle, is to be expected, due to the increase in the frequency of collisions that is typical of increasing pressures. However, more simulations covering a wider pressure range would be needed to confirm and reproduce this behavior. For the range where the permeances are weakly dependent on pressure, we calculated the average of the gas permeance for each gas and plotted it as a function of temperature ( Figure 3B). The plot shows that the average gas permeance decreases as the temperature increases for all the molecules. This behavior is likely to be a consequence of the increased kinetic energy of the molecules; by possessing higher kinetic energies (velocities), molecules more efficiently escape the attraction forces of the membrane. Therefore, the high temperature decreases the gas permeance by contrasting the attraction effects that steer the molecules toward the membrane. The low average permeance of N 2 in the entire range of temperatures considered in the simulations indicates that the attractive forces between the membrane and N 2 are considerably weaker, compared to those experienced by CO 2 and H 2 O. This fact is reflected in the potential energy profiles reported by Bartolomei and co-workers [27,49,50], describing the membrane-molecule interactions, showing that N 2 has the smallest well depth, followed by H 2 O and CO 2 . In spite of having the deepest potential energy well, CO 2 demonstrates permeance values lower than those of H 2 O. To explain such apparent illogicality requires more insight into the permeation process, schematized in Figure 4. However, upon closer scrutiny, the data in Figure 3A, although shown in the limited pressure range considered, seem to suggest a minimum level of permeance, followed by a rebound profile. Such a trend, in principle, is to be expected, due to the increase in the frequency of collisions that is typical of increasing pressures. However, more simulations covering a wider pressure range would be needed to confirm and reproduce this behavior. For the range where the permeances are weakly dependent on pressure, we calculated the average of the gas permeance for each gas and plotted it as a function of temperature ( Figure 3B). The plot shows that the average gas permeance decreases as the temperature increases for all the molecules. This behavior is likely to be a consequence of the increased kinetic energy of the molecules; by possessing higher kinetic energies (velocities), molecules more efficiently escape the attraction forces of the membrane. Therefore, the high temperature decreases the gas permeance by contrasting the attraction effects that steer the molecules toward the membrane. The low average permeance of N2 in the entire range of temperatures considered in the simulations indicates that the attractive forces between the membrane and N2 are considerably weaker, compared to those experienced by CO2 and H2O. This fact is reflected in the potential energy profiles reported by Bartolomei and co-workers [27,49,50], describing the membrane-molecule interactions, showing that N2 has the smallest well depth, followed by H2O and CO2. In spite of having the deepest potential energy well, CO2 demonstrates permeance values lower than those of H2O. To explain such apparent illogicality requires more insight into the permeation process, schematized in Figure 4. . An illustration of the permeation process of CO2, displayed through a pore unit of the graphtriyne layer (adopted from the potential energy curve reported in Ref. [50]). . An illustration of the permeation process of CO 2 , displayed through a pore unit of the graphtriyne layer (adopted from the potential energy curve reported in Ref. [50]).
First, the molecules located near the membrane are affected by the attraction forces and move toward the pores. In this step, the CO 2 molecules show a propensity toward more efficient adsorption by the membrane, due to its strong attraction forces. Then, a gas molecule can reach the center of the pore (z = 0), a stable configuration where the potential energy is minimal. If there is enough kinetic energy to overcome the attraction forces, the molecule will cross the membrane and fly to the other side. On the other hand, if the molecule does not have enough energy, it will remain adsorbed by the membrane until it acquires extra energy (collisions with other molecules) or assumes an optimal orientation for crossing. Due to the deeper potential well, CO 2 molecules need greater kinetic energy than H 2 O to overcome the attraction forces and to successfully cross the membrane. Moreover, in the case of CO 2 , the permeation is stereoselective with respect to the molecular orientation. The penetration process encounters energy barriers when CO 2 molecules approach the membrane in parallel configurations (see Figure 1 in Ref. [50] for the energy profiles corresponding to the perpendicular and parallel approach, and for the related discussions). This qualitatively explains why CO 2 molecules find it more difficult to pass through the membrane, even though the numbers of CO 2 molecules attracted and adsorbed by the membrane are larger than those of H 2 O (this can be seen in the z-density profiles, which will be discussed later). Therefore, as can be seen in Figure 3 for the single-layer systems, H 2 O permeances are the highest ones.
The gas permeances of the multilayer graphtriyne can be described using the data presented in Figure 5. We can see that the N 2 average permeances increase as the number of layers increases ( Figure 5B). This behavior is consistent with the increasingly stronger attractive forces probed by N 2 that the graphtryine is subjected to, passing from bi-to trilayers, as can be seen from the interaction energies of N 2 with the multilayers reported in Ref. [49]. However, CO 2 shows different behavior since the highest permeances, which are achieved in the bilayer system, decrease in the trilayer membrane ( Figure 5A). This is the consequence of the efficient adsorption of CO 2 in the interlayer regions between the three graphtriyne sheets (this aspect will be discussed later). The presence of strongly adsorbed molecules prevents other CO 2 molecules from crossing the membranes. Although the trend of permeances is consistent with the ones obtained in Ref. [50], the values reported here are lower (ranging from 0.6 to 3.0 × 10 7 GPU and from 0.3 to 0.6 × 10 7 GPU for CO2 and N2, respectively). Lower overall values of the gas permeances are a natural consequence of the presence of three different species comprising CO2/N2/H2O, competing to interact with the membrane in the gaseous mixtures (while the mixtures in Ref. [50] were binary) and by the higher temperatures characterizing the present simulations. In the case of H2O (Figure 5C), the average gas permeances do not change greatly according to the number of layers. Moreover, H2O permeances are significantly higher than those of CO2 and N2 for the single and trilayer systems (much less at low temperatures (see Figure 5)), with the noticeable exception of the bilayer system. Figure 5 also shows that the gas permeances for all the molecules decrease with increasing temperature for single-layer, bilayer, and trilayer systems (see also Figure 3B).
By comparison of the gas average permeances, we calculated the permeance selectivity values (obtained, for a pair of molecules, as the ratio of the individual gas permeances) and plotted them as a function of temperature for all types of membrane systems and pairs of molecules. The results are reported in Figure 6. It can be seen that the permeance selectivity is affected by temperature, so that we can assume such a Although the trend of permeances is consistent with the ones obtained in Ref. [50], the values reported here are lower (ranging from 0.6 to 3.0 × 10 7 GPU and from 0.3 to 0.6 × 10 7 GPU for CO 2 and N 2 , respectively). Lower overall values of the gas permeances are a natural consequence of the presence of three different species comprising CO 2 /N 2 /H 2 O, competing to interact with the membrane in the gaseous mixtures (while the mixtures in Ref. [50] were binary) and by the higher temperatures characterizing the present simulations. In the case of H 2 O (Figure 5C), the average gas permeances do not change greatly according to the number of layers. Moreover, H 2 O permeances are significantly higher than those of CO 2 and N 2 for the single and trilayer systems (much less at low temperatures (see Figure 5)), with the noticeable exception of the bilayer system. Figure 5 also shows that the gas permeances for all the molecules decrease with increasing temperature for single-layer, bilayer, and trilayer systems (see also Figure 3B).
By comparison of the gas average permeances, we calculated the permeance selectivity values (obtained, for a pair of molecules, as the ratio of the individual gas permeances) and plotted them as a function of temperature for all types of membrane systems and pairs of molecules. The results are reported in Figure 6. It can be seen that the permeance selectivity is affected by temperature, so that we can assume such a dependence, with the caveat that with the variations being lower than the standard error, the effect could be weak or poorly significant in the observed temperature range. The permeance selectivity of the H 2 O/N 2 pair is the highest one in both single-layer and trilayer systems, while it essentially overlaps that of CO 2 /N 2 in the bilayer system. Interestingly, H 2 O/CO 2 selectivity is~1 or slightly higher for both bi-and trilayer systems, indicating that such membranes are not selective for CO 2 -H 2 O separation. Figure 6 also shows that, in general, the trilayer system (C) has permeance selectivity values that are lower than for the single (A) and bilayer (B) systems. These are related to efficient molecular adsorption in the interlayer regions of the trilayer system, which can lower gas permeance (this aspect will be discussed in the next section).  Although the CO2/N2 permeance selectivity values obtained here are much lower than those reported for nanoporous graphene at 300 K by Liu et al. [23] (about 100, with CO2 permeance = 2.8 × 10 5 GPU) and by Schrier [17] (about 60, with CO2 permeance = 3 × 10 5 GPU) for porous graphene-E-stilbene-1 (PG-ES1) at 325 K; the CO2 permeances for bilayer graphtriyne (ranging from 0.6 to 3.0 × 10 7 GPU) are two orders of magnitude larger. Moreover, the CO2/N2 permeance selectivity values for the bilayer system are comparable to and higher than the value (5.4) reported in our previous work for CO2/N2 binary mixtures at 300 K [50]. The CO2/N2 permeance selectivity values are also comparable with those reported by Wu and co-workers [24] for fluorine-modified nanoporous graphene at 300 K (ranging from 4 to 11). This indicates that the bilayer graphtriyne membrane could represent an efficient molecular sieve to be used in an initial separation step for CO2 postcombustion separation; the flue gas is mainly composed of molecular nitrogen.

Gas Adsorption
The adsorption process can be characterized by the mean number densities of gas molecules along the -axis (the direction perpendicular to the graphtriyne membrane). Although the CO 2 /N 2 permeance selectivity values obtained here are much lower than those reported for nanoporous graphene at 300 K by Liu et al. [23] (about 100, with CO 2 permeance = 2.8 × 10 5 GPU) and by Schrier [17] (about 60, with CO 2 permeance = 3 × 10 5 GPU) for porous graphene-E-stilbene-1 (PG-ES1) at 325 K; the CO 2 permeances for bilayer graphtriyne (ranging from 0.6 to 3.0 × 10 7 GPU) are two orders of magnitude larger. Moreover, the CO 2 /N 2 permeance selectivity values for the bilayer system are comparable to and higher than the value (5.4) reported in our previous work for CO 2 /N 2 binary mixtures at 300 K [50]. The CO 2 /N 2 permeance selectivity values are also comparable with those reported by Wu and co-workers [24] for fluorine-modified nanoporous graphene at 300 K (ranging from 4 to 11). This indicates that the bilayer graphtriyne membrane could represent an efficient molecular sieve to be used in an initial separation step for CO 2 post-combustion separation; the flue gas is mainly composed of molecular nitrogen.

Gas Adsorption
The adsorption process can be characterized by the mean number densities of gas molecules along the z-axis (the direction perpendicular to the graphtriyne membrane). For the single-layer system, these z-density profiles, ρ(z), tend to peak around a distance of 3.4 Å from the surface in all gas molecules ( Figure 7A). The strength of the attraction exerted by the membrane affects the height of the z-density profile peaks; accordingly, CO 2 shows the highest peaks. We can also see from the radial distribution functions g(r) ( Figure 7B) that CO 2 has the highest probability to be found near the carbon atoms of the graphtriyne (around 4 Å). The strong peaks of CO 2 in the z-density and the radial distribution function are a consequence of the deep and wide potential well that characterizes the CO 2 -membrane interaction potential profile, resulting in a strong long-range attraction (see Figure 4 and Refs. [49,50]). As already discussed in the previous section, the permeation events are closely related to the adsorption of gas molecules over the surfaces of the membranes. The more the molecules are adsorbed, the higher their permeance. However, the stereodynamic requirements of CO 2 for the membrane-crossing process (which is reported and exhaustively discussed in our previous work [50]) lead to permeances lower than those of H 2 O.
Molecules 2022, 27, x FOR PEER REVIEW 9 of 16 the surfaces of the membranes. The more the molecules are adsorbed, the higher their permeance. However, the stereodynamic requirements of CO2 for the membrane-crossing process (which is reported and exhaustively discussed in our previous work [50]) lead to permeances lower than those of H2O. The gas uptake can be estimated by integrating the area under the peaks of thedensity profiles. The interval of integration along the z-direction is the adsorption region located within 6.9 Å with respect to the membrane position of the single-layer system. Notice that the -density does not peak at the center of the pore = 0 , where, according to the potential profiles (see Figure 4), we see the energy minimum. This feature ofdensity has already been discussed in our previous report [50] and has been interpreted, in terms of the oscillatory motion of the molecule in a physisorption state, as being due to the longer residence times at turning points and, in the case of CO2, are also enhanced by the stereodynamics of the crossing process [49].
The gas uptakes for all applied pressures are reported in Figure 8A in the form of adsorption isotherms. We can see that the gas uptake is linearly related to the initial pressures and increases with the increasing values of pressure. Therefore, using linear regression, we can estimate the adsorption coefficients as the slopes of the linear functions fitting the adsorption isotherms. The adsorption coefficients for all gas molecules in the single-layer system are reported in Figure 8B, as a function of the temperature. It can be seen that CO2 has the highest gas uptakes and adsorption coefficients among the gas molecules. On the other hand, despite the fact that H2O has the highest permeance, its gas uptakes and adsorption coefficients are lower, again indicating weaker attraction forces than those acting for CO2. The weak attraction of N2 by the membrane is manifested by the fact that N2 has the lowest adsorption coefficient (about 0.06 to 0.07 mmol g −1 atm −1 ). Figure 8B shows that the adsorption coefficients decrease for increasing temperature. This trend appears to fluctuate at 350 K; this still leaves room for a different trend than the expected linear one. However, a much longer simulation and a wider temperature range are needed to obtain a clear view. Obviously, with higher temperatures, the gas molecules The gas uptake can be estimated by integrating the area under the peaks of the zdensity profiles. The interval of integration along the z-direction is the adsorption region located within ± 6.9 Å with respect to the membrane position of the single-layer system. Notice that the z-density does not peak at the center of the pore (z = 0), where, according to the potential profiles (see Figure 4), we see the energy minimum. This feature of z-density has already been discussed in our previous report [50] and has been interpreted, in terms of the oscillatory motion of the molecule in a physisorption state, as being due to the longer residence times at turning points and, in the case of CO 2 , are also enhanced by the stereodynamics of the crossing process [49].
The gas uptakes for all applied pressures are reported in Figure 8A in the form of adsorption isotherms. We can see that the gas uptake is linearly related to the initial pressures and increases with the increasing values of pressure. Therefore, using linear regression, we can estimate the adsorption coefficients as the slopes of the linear functions fitting the adsorption isotherms. The adsorption coefficients for all gas molecules in the single-layer system are reported in Figure 8B, as a function of the temperature. It can be seen that CO 2 has the highest gas uptakes and adsorption coefficients among the gas molecules. On the other hand, despite the fact that H 2 O has the highest permeance, its gas uptakes and adsorption coefficients are lower, again indicating weaker attraction forces than those acting for CO 2 . The weak attraction of N 2 by the membrane is manifested by the fact that N 2 has the lowest adsorption coefficient (about 0.06 to 0.07 mmol g −1 atm −1 ). Figure 8B shows that the adsorption coefficients decrease for increasing temperature. This trend appears to fluctuate at 350 K; this still leaves room for a different trend than the expected linear one. However, a much longer simulation and a wider temperature range are needed to obtain a clear view. Obviously, with higher temperatures, the gas molecules have a higher tendency to escape; thus, the physisorption is expected to be less effective. Similar trends are observed in the case of bilayer and trilayer systems. CO2 shows the highest gas uptake and adsorption coefficients ( Figure 9B,D), and the adsorption coefficients also depend on the temperature: higher temperatures provide lower adsorption coefficients. Unlike the single-layer system, in the bilayer and trilayer systems, we find -density peaks in the interlayer regions ( Figure 9A,C). Therefore, it is expected that the bilayer and trilayer membranes should be more efficient in adsorbing molecules due to their interlayer pores (as reported in the Supplementary Information, Tables S11-S14). However, adsorbed molecules can, in some conditions, also contribute to diminishing the gas permeance by saturating the pores and preventing other molecules from crossing the membrane. This phenomenon can be observed in the case of CO2 with trilayer systems, due to the strong attraction and high uptake. Indeed, we can see from Figure 5A that the CO2 permeance values significantly increase when passing from singlelayer to bilayer systems, then collapsing to values typical of single layers when passing from bilayer to trilayer systems. This phenomenon can then accumulate into a decrease in the value of the permeance selectivity of the trilayer system, relative to other systems, as shown in Figure 6C. Similar trends are observed in the case of bilayer and trilayer systems. CO 2 shows the highest gas uptake and adsorption coefficients ( Figure 9B,D), and the adsorption coefficients also depend on the temperature: higher temperatures provide lower adsorption coefficients. Unlike the single-layer system, in the bilayer and trilayer systems, we find z-density peaks in the interlayer regions ( Figure 9A,C). Therefore, it is expected that the bilayer and trilayer membranes should be more efficient in adsorbing molecules due to their interlayer pores (as reported in the Supplementary Information, Tables S11-S14). However, adsorbed molecules can, in some conditions, also contribute to diminishing the gas permeance by saturating the pores and preventing other molecules from crossing the membrane. This phenomenon can be observed in the case of CO 2 with trilayer systems, due to the strong attraction and high uptake. Indeed, we can see from Figure 5A that the CO 2 permeance values significantly increase when passing from single-layer to bilayer systems, then collapsing to values typical of single layers when passing from bilayer to trilayer systems. This phenomenon can then accumulate into a decrease in the value of the permeance selectivity of the trilayer system, relative to other systems, as shown in Figure 6C. The interlayer region is relatively selective to CO2 molecules (see Figure 9A,C); for instance, we obtained interlayer adsorption selectivity of CO2/N2 = 20.23, CO2/H2O = 1.85 for the bilayer system at 333 K and 4.00 atm, with CO2/N2 = 42.45, CO2/H2O = 2.38 for the trilayer system at 400 K and 2.54 atm (see the Supplementary Information, Table S18). The interlayer adsorption refers to the adsorption values in the interlayer region, i.e., between the two outermost red lines of the z-density profiles (see Figure 9A,C). The high interlayer selectivity is an additional feature of multilayer graphtriyne membranes. However, for most practical applications, it is better to express the selectivity in terms of the total adsorption selectivity. This is because the available experimental and theoretical data reported in the literature often only report the total adsorptions (in pores + surface adsorptions). The total adsorption refers to the sum of the adsorption values in the interlayer region and in the outer (adsorption) region, within 6.9 Å from the membrane position. The interlayer and total adsorption selectivity relative to two molecules, A and B.
/ , was calculated according to the formula: where and are the numbers of adsorbed molecules A and B, respectively, while and are the numbers of free molecules of A and B, respectively [63]. The total adsorption selectivity values for all the systems are reported in Table 1, while some of the data are plotted in Figure 10A. Figure 10 and Table 1 show that trilayer graphtriyne exhibits relatively high CO2/N2 total adsorption selectivity, especially at low pressure and temperature. As an example, at a temperature of 353 K, the CO2/N2 total adsorption selectivity ranges from 4.10 to 19.45, while the CO2/H2O total adsorption selectivity ranges from 1.59 to 11.08, depending on the applied initial pressure. The interlayer region is relatively selective to CO 2 molecules (see Figure 9A,C); for instance, we obtained interlayer adsorption selectivity of CO 2 /N 2 = 20.23, CO 2 /H 2 O = 1.85 for the bilayer system at 333 K and 4.00 atm, with CO 2 /N 2 = 42.45, CO 2 /H 2 O = 2.38 for the trilayer system at 400 K and 2.54 atm (see the Supplementary Information, Table S18). The interlayer adsorption refers to the adsorption values in the interlayer region, i.e., between the two outermost red lines of the z-density profiles (see Figure 9A,C). The high interlayer selectivity is an additional feature of multilayer graphtriyne membranes. However, for most practical applications, it is better to express the selectivity in terms of the total adsorption selectivity. This is because the available experimental and theoretical data reported in the literature often only report the total adsorptions (in pores + surface adsorptions). The total adsorption refers to the sum of the adsorption values in the interlayer region and in the outer (adsorption) region, within ±6.9 Å from the membrane position. The interlayer and total adsorption selectivity relative to two molecules, A and B. S A/B ads , was calculated according to the formula: where n A (ads) and n B (ads) are the numbers of adsorbed molecules A and B, respectively, while n A ( f ree) and n B ( f ree) are the numbers of free molecules of A and B, respectively [63]. The total adsorption selectivity values for all the systems are reported in Table 1, while some of the data are plotted in Figure 10A. Figure 10 and Table 1 show that trilayer graphtriyne exhibits relatively high CO 2 /N 2 total adsorption selectivity, especially at low pressure and temperature. As an example, at a temperature of 353 K, the CO 2 /N 2 total adsorption selectivity ranges from 4.10 to 19.45, while the CO 2 /H 2 O total adsorption selectivity ranges from 1.59 to 11.08, depending on the applied initial pressure. These CO 2 /N 2 total adsorption selectivity values are comparable to those of some porous materials for CO 2 capture applications, such as covalent organic polymers (8.4-13.7 at 1.01 bar, 298 K) [18], metal-organic frameworks (5.0-40.0 at 20 bar, 298 K) [14], and functionalized graphitic slit pores (5.0-20.0 at 20 bar, 298 K) [63]. We can also visually verify the selectivity by looking at the snapshot of a configuration sampled from a typical simulation run, as shown in Figure 10B, where the CO 2 molecules that are adsorbed are predominant. In general, the total adsorption selectivity of all the gas pairs shows a decreasing trend as the pressure and temperature increase (Table 1 and Figure 10A).  These CO2/N2 total adsorption selectivity values are comparable to those of some porous materials for CO2 capture applications, such as covalent organic polymers (8.4-13.7 at 1.01 bar, 298 K) [18], metal-organic frameworks (5.0-40.0 at 20 bar, 298 K) [14], and functionalized graphitic slit pores (5.0-20.0 at 20 bar, 298 K) [63]. We can also visually verify the selectivity by looking at the snapshot of a configuration sampled from a typical simulation run, as shown in Figure 10B, where the CO2 molecules that are adsorbed are predominant. In general, the total adsorption selectivity of all the gas pairs shows a decreasing trend as the pressure and temperature increase (Table 1 and Figure 10A).
As a matter of fact, the CO2/N2 total selectivity values reported in Table 1 are higher than those given in our previous work [50]. For instance, here, we obtained the values of the CO2/N2 total adsorption selectivity for the trilayer system, ranging from 4.8 to 8.9 (trinary gaseous mixture at 333 K), while in our previous work, we reported values ranging from 4.0 to 8.0 (binary gaseous mixture at 300 K). However, the gas uptake values and the corresponding adsorption coefficients are relatively lower than those obtained in the previous work. It is obvious that this difference comes from the different kinds of gaseous mixtures and conditions applied in the simulations. Nevertheless, both reports are in agreement, predicting that the trilayer membrane exhibits a high selectivity for CO2 capture. Therefore, when presenting all the advantages of carbon-based materials, multilayer graphtriyne membranes are promising alternatives for post-combustion CO2 capture and separation. In particular, the bilayer graphtriyne membrane, with its permeability and permeance selectivity, exhibits good performance as an initial molecular sieve candidate for post-combustion CO2 separation, whereas trilayer graphtriyne, with its high gas uptake and adsorption selectivity, is comparable and competitive with other carbon-based adsorbing materials for post-combustion CO2 capture.

Conclusions
In this work, we performed extensive molecular dynamics simulations to assess the suitability of multi-layer graphtriyne membranes for CO2 capture and separation. To improve the accuracy of the results, we adopted a proper formulation for the involved force fields, which represent, in an internally consistent way, the fundamental interaction components for both molecule-molecule and molecule-layer (surface) systems. Comprehensive sets of values for the permeance, uptake, and selectivity coefficients have been obtained, which also provide a dynamic and stereodynamic interpretation of the observed data in an interval of temperatures and pressures of interest for many applications. Remarkably, this data set can justifiably be considered predictive in terms of the care taken in the description of the atomistic molecular dynamics and intermolecular interactions, as derived in previous works. Having uniformly distributed and tunable As a matter of fact, the CO 2 /N 2 total selectivity values reported in Table 1 are higher than those given in our previous work [50]. For instance, here, we obtained the values of the CO 2 /N 2 total adsorption selectivity for the trilayer system, ranging from 4.8 to 8.9 (trinary gaseous mixture at 333 K), while in our previous work, we reported values ranging from 4.0 to 8.0 (binary gaseous mixture at 300 K). However, the gas uptake values and the corresponding adsorption coefficients are relatively lower than those obtained in the previous work. It is obvious that this difference comes from the different kinds of gaseous mixtures and conditions applied in the simulations. Nevertheless, both reports are in agreement, predicting that the trilayer membrane exhibits a high selectivity for CO 2 capture. Therefore, when presenting all the advantages of carbon-based materials, multilayer graphtriyne membranes are promising alternatives for post-combustion CO 2 capture and separation. In particular, the bilayer graphtriyne membrane, with its permeability and permeance selectivity, exhibits good performance as an initial molecular sieve candidate for post-combustion CO 2 separation, whereas trilayer graphtriyne, with its high gas uptake and adsorption selectivity, is comparable and competitive with other carbon-based adsorbing materials for post-combustion CO 2 capture.

Conclusions
In this work, we performed extensive molecular dynamics simulations to assess the suitability of multi-layer graphtriyne membranes for CO 2 capture and separation. To improve the accuracy of the results, we adopted a proper formulation for the involved force fields, which represent, in an internally consistent way, the fundamental interaction components for both molecule-molecule and molecule-layer (surface) systems. Comprehensive sets of values for the permeance, uptake, and selectivity coefficients have been obtained, which also provide a dynamic and stereodynamic interpretation of the observed data in an interval of temperatures and pressures of interest for many applications. Remarkably, this data set can justifiably be considered predictive in terms of the care taken in the description of the atomistic molecular dynamics and intermolecular interactions, as derived in previous works. Having uniformly distributed and tunable pores, with all the advantages of carbon-based materials, multilayer graphtriyne membranes are promising candidates for the separation and storage of CO 2 from post-combustion flue gases composed of CO 2 , N 2 and H 2 O gaseous mixtures. The bilayer graphtriyne membrane represents a good alternative as a molecular sieve for CO 2 separation, while the trilayer graphtriyne membrane is promising for post-combustion CO 2 storage and is competitive compared to other carbon-based adsorbing materials. Further developments include the extension of the approach to new carbon-based materials and molecules and the construction of a general force field for accurate simulations of gas-membrane systems.