In Silico Design and Characterization of Graphene Oxide Membranes with Variable Water Content and Flake Oxygen Content.

Graphene oxide (GO) membranes offer exceptional promise for certain aqueous separation challenges, such as desalination. Central to unlocking this promise and optimizing performance for a given separation is the establishment of a detailed molecular-level understanding of how the membrane's composition affects its structural and transport properties. This understanding is currently lacking, in part due to the fact that, until recently, molecular models with a realistic distribution of oxygen functionalities and interlayer flake structure were unavailable. To understand the effect of composition on the properties of GO membranes, models with water contents and oxygen contents, varying between 0% and 40% by weight, were prepared in this work using classical molecular dynamics simulations. The change in membrane interlayer distance distribution, water connectivity, and water diffusivity with water and oxygen content was quantified. Interlayer distance distribution analysis showed that the swelling of GO membranes could be controlled by separately tuning both the flake oxygen content and the membrane water content. Water-molecule cluster analysis showed that a continuous and fully connected network of water nanopores is not formed until the water content reaches ∼20%. The diffusivity of water in the membrane was also found to strongly depend on both the water and the oxygen content. These insights help understand the structure and transport properties of GO membranes with sub-nanometer interlayer distances and could be exploited to enhance the performance of GO membranes for aqueous separation applications. More broadly, the high-throughput in silico approach adopted could be applied to other nanomaterials with intrinsic non-stoichiometry and structural heterogeneity.

Graphene oxide (GO) flakes can be fabricated into multilayered and semipermeable materials known as GO membranes.[1] The individual nanopores of GO membranes are highly permeable to water but not necessarily dissolved solutes,[2,3] and this is a particularly desirable characteristic for an aqueous separation membrane.[4] To improve the efficacy of a given membrane separation, water permeability, and water-solute selectivity should both be optimized, but these two properties are typically negatively correlated. For example, in the case of desalination, reducing the average interlayer distance of GO membranes improves the rejection rate of dissolved salts but decreases the water permeability.[5,6] Central to the challenge of optimizing the separation performance of GO membranes is improving the understanding of the complex interplay between a membrane’s chemical composition, interlayer structure, and nanopore morphology.

An important feature of GO membranes is that water readily adsorbs into the nanopores due to the inherent hydrophilicity of the constituent GO flakes. The interlayer distance is ∼0.7 nm in dry conditions but increases with relative humidity, and when membranes are immersed in liquid water, the interlayer distance swells to ∼1.3 nm.[7−9] Tuning the membrane water content therefore provides a convenient means of controlling the interlayer distance and the effective nanopore diameter. Several non-destructive methods have been investigated for the control of water content in GO membranes, such as physical confinement,[6] solvent drying[10,11] or cross-linking with multivalent ions.[12,13] A significant improvement in salt rejection has been achieved by confining membranes within an epoxy resin, which serves to restrict the swelling.[6] This improvement is due to an increase in the ion-permeation free-energy barrier, associated with an enhancement of ion dehydration effects due to increasing nanoscale confinement upon the ion entering the membrane nanopores.[14]

The humidity-dependent interlayer distance of GO membranes differs significantly between samples, and this is largely attributable to variability in the chemical composition of the constituent GO flakes.[7] The most common oxidizing methods used to prepare GO flakes result in different total oxygen content and types of oxygen-containing functionalities.[15,16] For example, the Brodie oxidation method[17] is known to result in a lower overall oxygen content than the more commonly employed Hummers method[18] and favors the formation of conjugated epoxide and hydroxide functionalities.[15] As a result, differences in water permeability have been observed[19,20] depending on whether the flakes were oxidized using the Hummers,[18] Brodie,[17] or Staudenmaier[21] methods. The mechanical properties of GO membranes are also affected by changes in the degree and type of oxidation.[1,22,23] The total oxygen content of GO can also be controlled by chemical, thermal, and electrochemical reduction.[5,24,25] Reduced graphene oxide (RGO) membranes are chemically distinct and more hydrophobic[24,26] than their GO analogues, resulting in different transport properties,[19,20] including the possibility of improved salt rejection.[5,27]

One task of paramount importance is the systematic and quantitative characterization of GO membranes across the wide range of accessible membrane water contents and flake oxygen contents. This assignment is ideally suited to molecular simulation, which can be used for the rapid screening of candidate materials and characterization of properties that may be difficult to probe experimentally. In this study, the changes in GO membrane properties (interlayer distance, water connectivity, and water diffusivity) due to membrane composition (water content and flake oxygen content) were systematically investigated in silico using atomistic GO membrane models and classical MD simulations, with the ultimate aim of underpinning the design of GO membranes with improved performance for aqueous separations.

Results and Discussion

GO membrane models were constructed using the multistep molecular dynamics (MD) routine developed in our previous work.[28] In their review on the computer modeling of carbon-based porous materials, Striolo et al.(29) highlighted model structures of GO membranes as a key area in need of improvement to better understand their properties. In particular, they noted that often-ignored nanoscale features such as variations in pore size and the density and distribution of oxygen-containing groups are critical in understanding and predicting the permeability and selectivity of this class of materials. Our approach to constructing GO membrane models, outlined in the Methodology section, addresses this need. It also has a number of distinct advantages over other models described in the literature that render it suitable for the in silico design of GO membranes:

The GO flake model is fully atomistic and has oxygen-containing groups (epoxide, hydroxide, and carboxylic acid). These functionalities randomly decorate the plane and edges of a flake, ensuring oxygen-functionality distributions are not correlated, capturing the inhomogeneous and non-stoichiometric nature of oxidation.[15,30] An absence of oxygen functionalities in the GO nanopores will result in a spurious determination of the effective nanopore diameter and is known to result in significant differences in water dynamics.[31]

GO flakes are fully flexible. Qualitative differences have been observed in the phase behavior, structure, and diffusivity of nanopore water in graphene channels when modeled with rigid compared with flexible walls.[32]

Membrane models with a water content of 15% by weight and oxygen content of 25% by weight have an average interlayer distance, bulk density, pore volume, and elastic modulus in excellent agreement with experimental studies.[28]

The multilayered structure is generated by pulling and depositing GO flakes and water onto a surface from a random initial configuration, mimicking the real fabrication process (i.e., flow-directed self-assembly onto a filter or support).[1,3,33,34] The process distinguishes the resulting models from the chemically similar but structurally distinct graphite oxides, which inherit their long-range order from the parent graphite material.

The interlayer distance is not defined a priori, unlike in numerous other simulation studies of GO membranes,[6,14,31,35−41] but rather adopts a value, depending on membrane composition, after preparation and equilibration of the model.

The method is computationally inexpensive (∼1000–1500 CPU hours per model, depending on the composition). Although the computational characterization of crystalline nanoporous materials with well-defined structures and atomic coordinates (e.g., metal–organic frameworks[42] or zeolites)[43] can be relatively straightforward, it is significantly more challenging for materials that have intrinsic non-stoichiometry and structural heterogeneity, such as GO membranes.[8,15] One challenge is that the unambiguous specification of atomic coordinates in a GO flake is not possible. In addition, the interlayer distance and related properties (e.g., pore size distribution) may vary considerably throughout the membrane structure and cannot simply be calculated for a unit cell and extrapolated to experimentally relevant length scales, as is routinely performed for crystalline materials. The approach is therefore amenable to generating an ensemble of models from which the average properties and their uncertainty can be determined. The low computational cost of this approach is principally associated with simulating GO flake assembly in the presence of a small number of water molecules to exfoliate GO flakes, avoiding costly simulations of liquid water. In addition, it provides a convenient means to incorporate water into the membrane nanopores without having to perform computationally expensive Monte Carlo simulations, in which the membrane is in equilibrium with a water reservoir at a fixed relative humidity.

The membrane water content, mw, and flake oxygen content, mo, of membrane models are expressed as weight percentages:where mwater is the total mass of water molecules in the membrane, nGO is the number of constituent GO flakes in the membrane (nGO = 6 in this work), mGO is the mass of a GO flake, and moxygen is the total mass of oxygen on the flake present on the basal plane only. The mass of oxygen atoms on the edges of the flake are not included in the calculation of moxygen because their inclusion results in a non-representative low surface density of oxygen on the plane at low oxygen content. This is a finite size effect due to the increased relative contribution of the edges compared with real GO flakes, which have a much larger size. To separately investigate the effects of mw and mo on membrane properties these values were varied from 0% to 40% with 5% intervals. The range of water contents studied was chosen to reflect membranes with different degrees of hydration, ranging from the anhydrous state to the fully swollen state (i.e., when immersed in liquid water). For each membrane composition (mw and mo), 10 models were constructed, from which mean properties and uncertainties were determined. Therefore, a total of 180 GO membrane models were constructed in this study. Models with variable m were prepared with a constant oxygen content (m = 25%, C-to-O atomic ratio = 4:1). The range of oxygen contents was chosen to reflect the typical values of GO using standard oxidation methods (mo = 25% to 40%, C-to-O atomic ratio between 4:1 and 2:1)[15] and those accessible by reduction to RGO (mo < 25%).[26] Models with variable mo values were prepared with a constant number of water molecules (1067), which corresponds to mw = 15% for a membrane with mo = 25%, roughly typical for the water content of a membrane at ambient humidity.[44,45]

The interlayer distance of GO membranes is perhaps the most important property to characterize and quantify, because its value, minus the effective flake thickness, defines the approximate diameter of the nanopores through which mass transport occurs. Experimentally, swelling of GO membranes is most commonly investigated by analysis of X-ray diffraction (XRD) patterns, where the position of the (001) reflection is assumed to be equal to the average interlayer distance, dGO. In this work, dGO was calculated from dGO = L/nGO, where L is the final length of the simulation cell containing the membrane in the z dimension. Figure shows the change in dGO with water and oxygen content.

Figure 1

Average interlayer distance, dGO, for membranes prepared with variable membrane water content, mw, (black, circles) and flake oxygen content, m0 (red, squares). For most data points, the error bars are smaller than the symbols. The blue labels correspond to separate regions of the swelling profile for the variable water content models.[44]

For the dry membrane (mw = 0%), dGO = 0.77 nm, but the membrane swells to dGO = 1.35 nm at the highest water content studied (mw = 40%). Experimentally, the swelling of GO membranes is studied by increasing the relative humidity (RH), and two distinct regions of the swelling profile can typically be identified. In the first region (RH = 30–75%), the degree of swelling is modest, corresponding to an increase in dGO of ∼0.1 nm.[7,46−48] In the second region, at higher humidities (RH > 75%) or when immersed in liquid water,[46] a more significant increase in dGO of ∼0.3 nm is observed. Although the direct relationship between RH and dGO was not made in the present study, Korobov et al.(44) reported that a membrane with a RH of 30% and 75% had water contents of 17% and 24%, respectively. This correspondence between RH and water content allows us to approximately separate the simulated swelling profile into three regions for the convenience of the following discussion; regions A, B, and C, defined by m = 0% to 17%, 17% to 24%, and 24% to 40%, respectively. In region A, which corresponds to RH < 30%, the average interlayer distance is approximately constant (dGO ≈ 0.80 nm). In region B, which corresponds to 30% < RH < 75%, there is a modest increase in average interlayer distance (dGO = 0.80 to 0.96 nm). In region C, which corresponds to RH > 75%, the average interlayer distance increases more rapidly (dGO = 0.96 to 1.35 nm), in qualitative agreement with the experimental literature. It is interesting to compare these observations with the simulated results of Devanathan et al.,[45] which showed an abrupt increase in the interlayer distance, from 0.8 nm at mw = 0% to 1.1 nm at mw = 1% and then a constant interlayer distance for all other water contents (up to mw = 23%). This discontinuous increase in interlayer distance at low water content is inconsistent with our results and the available experimental evidence,[5,7,46−49] which show a continuous and gradual expansion of interlayer distance upon hydration.

Region A is not reported in experiments, and this is probably related to the known difficulty in removing all of the water from the membrane;[48] thus, membrane models with very low water content may not have an easily realizable experimental counterpart. The changes in interlayer distance with water content can be explained in the context of the membranes’ bulk density. The density passes through a maximum at a water content of 15% (Figure S2), and this is due to the inherent porosity of the membrane. At low water content (region A), water can be added into the free pore volume without a change in the interlayer distance, increasing the bulk density. At higher water contents (regions B and C), the bulk density decreases. This is because more water can only be accommodated into the pores by increasing the interlayer distance (swelling) and because the skeletal density of the GO framework (∼2.2 g cm–3) is much higher than the bulk density of liquid water (1.0 g cm–3).

Figure shows that the average interlayer distance of membranes with variable oxygen content increases linearly from dGO = 0.68 to 0.90 nm because an increase in the number of oxygen functionalities increases the effective flake thickness. This suggests that dGO could be controlled by tuning the chemical composition of the constituent GO flakes by oxidation and reduction reactions. For membranes with mo = 0%, dGO is much greater than the carbon–carbon distance between adjacent layers in graphite (0.34 nm), and this difference is due to the presence of water and imperfect packing of flakes in the model. Experimentally, a variation in average interlayer distance has been reported for GO samples with different degrees and types of oxygen functionality[7] and control of the interlayer distance by tuning oxygen content has already been realized.[5,50] Cho et al.(5) varied the oxygen content of vacuum-dried membranes using different oxidation methods and thermal treatment. Their XRD analysis showed that Hummers membranes, with a C-to-O ratio of 2.0 (mo = 40.0%), had an interlayer distance of 0.85 nm, and Brodie membranes, with a C-to-O ratio of 2.1 (m = 38.9%), had a slightly lower interlayer distance of 0.82 nm. The surface of GO flakes oxidized using the Hummers or Brodie method is known to be dominated by hydroxide and epoxide functionalities[23] as well as some isolated regions of aromatic carbon,[51] consistent with the structure of our flake models (Figure S1).[28] The value of 0.85 nm is less than that predicted by our models at the same oxygen content, but an unambiguous quantitative comparison to the experiment is not possible because the exact water content of these membranes is unknown.

Further insight into the change in interlayer structure of GO membranes can be obtained by visualizing individual models. One set of models with increasing membrane water content is shown in Figure . From Figure , it is clear that GO membranes have considerable variation in interlayer distance between different regions of the membrane as well as between a given pair of adjacent GO flakes. The layering of GO flakes is far from ideal, which is not unexpected given the broad XRD peaks observed experimentally.[8,52,53] This heterogeneity means that simple models of GO nanopores, defined by perfectly parallel and graphene plates with a constant interlayer distance are structurally unrepresentative. The variation in interlayer distance revealed by Figure suggests it would be more instructive to study the distribution of interlayer distances, P(d), rather than the averaged value, dGO, as in Figure . Interlayer distance distributions for membranes with different water content and flake oxygen content are shown in Figure .

Figure 2

Swelling of GO membranes with increasing water content. GO flake carbon, oxygen, and hydrogen atoms are shown with black, red, and yellow spheres, respectively; the cyan surface indicates a molecular surface generated over nanopore water molecules.

Figure 3

Distribution of interlayer distances, P(d), for GO membranes with variable (a) water content, mw, and (b) oxygen content, mo.

Figure a shows how the distribution of interlayer distances changes with water content. The dry membrane has three well-defined peaks, at d = 0.49, 1.00, and 1.51 nm, and these peaks indicate strong correlations of GO flake pairs in direct contact, separated by one other GO flake and separated by two other GO flakes, respectively. The position of the first peak is much less than the average interlayer distance for the same membrane (dGO = 0.77 nm), which is due to the inefficient packing of flakes and the presence of voids in the structure when prepared without water (as shown in Figure ). The correlation between flake positions decreases with distance leading to a reduction in the intensity of successive P(d) peaks.

As the water content is increased, the intensity of the first peak decreases and shifts to slightly larger distances due to the incorporation of water molecules into the interlayer space. At mw = 15%, a new peak becomes evident (at d = 0.72 nm), corresponding to an interlayer distance where adjacent GO flakes are separated by a water monolayer. The peaks at 0.49 and 0.72 nm are similar in intensity at this water content, demonstrating coexistence of dehydrated regions and monolayer water regions. At higher water contents, the peak at 0.72 nm becomes more intense as the number and size of monolayer water regions increases and the peak initially at 0.49 nm eventually completely disappears (mw > 30%). At these higher water contents, all GO flakes are separated by at least one water layer, and regions in which flakes are in direct contact with each other are no longer observed. At the highest water content, three well-defined peaks in P(d) can be identified (at d = 0.73, 1.38, and 2.11 nm), corresponding to regions in which GO flakes are separated by a monolayer, bilayer, and trilayer of water, respectively. The value of dGO for this membrane is 1.35 nm (Figure ), which is very close to the position of the bilayer peak. Despite this, many monolayer and trilayer water regions are also present at this water content, and this heterogeneity will have important implications for understanding the permeability of GO membranes in the very swollen state.

It is helpful to validate these swelling observations by discussion of the available experimental evidence. Lerf et al.(7) investigated the in situ hydration of graphite oxide using XRD analysis. They observed a single (001) peak in the XRD pattern for an anhydrous sample. Within the first hour of exposure to 100% RH conditions, they observed the gradual disappearance of this peak and appearance of a second peak in the XRD pattern, and these two peaks coexist at certain intermediate times. The coexistence of two distinct interlayer distances upon hydration can be attributed to interstratification (i.e., randomly packed hydration states throughout the interlayer space of the membrane).[46] Finally, at longer time scales, they observed a continuous shift of the second reflection in the XRD pattern to larger distances, as well as a broadening of the peak, corresponding to an increase in the number of different interlayer distances. Scanning from left to right in Figure a can therefore be interpreted as equivalent to the time-dependent swelling of GO membranes in the study by Lerf.[7] Similar observations were made in a separate study.[48] In support of this finding, scanning force microscopy imaging of bilayer GO showed that, up to RH = 80%, the interlayer space gradually increases by incorporation of immobile patches of water less than 10 nm in size[46] and not by a sudden, discontinuous incorporation of well-defined integer layers of water throughout the entire membrane. Furthermore, AFM analysis[54] of a GO membrane sample with an average interlayer distance 0.78 nm showed a trimodal interlayer distance distribution, which was reproduced by fitting to 3 Gaussian functions, centered at 0.48, 0.76, and 1.11 nm. This compares well with the peak positions in P(d) for the mw = 15% model, which has a similar average interlayer distance (0.80 nm), and P(d) peaks positioned at 0.52, 0.72, and 1.08 nm.

Figure b shows the change in P(d) for models with increasing oxygen content. For mo = 0%, the distribution consists of a series of well-defined intense peaks, the first 3 of which are positioned at 0.37, 0.75, and 1.14 nm. The position of the first peak and distance between this and subsequent peaks are close to the interlayer distance of graphite because there are no oxygen functionalities on the plane of the flake in this model. The absence of a peak associated with a water monolayer implies that the water does not reside in the interlayer space in these membranes. Water instead accumulates at the hydrophilic flake edges, which are functionalized with carboxylic acid groups rather than the hydrophobic planes. This effect can be visualized in Figure and is supported by an analysis of water molecule coordination numbers (Figure S6), which show a decrease in the water–water and water–carboxylic acid coordination numbers with oxygen content, associated with water being increasingly present in the interlayer space, instead of at the flake edges. This effect also explains the presence of some larger nanopores shown in the pore size distributions for membranes with very low oxygen content (Figure S4).

Figure 4

GO membranes with variable flake oxygen contents. GO flake carbon, oxygen, and hydrogen atoms are shown in black, red, and yellow, respectively; the cyan surface indicates a molecular surface generated from nanopore water.

As mo is increased, the intensity of the first peak decreases and gradually shifts to larger distances. At intermediate oxygen contents the structure is poorly defined beyond the first peak, indicating poor correlations between GO flake interlayer distances. At high oxygen content (mo = 40%), the peak positions again become well-defined, with distinct peaks at 0.63, 1.23, and 1.80 nm. This can be explained in terms of flake surface roughness. At intermediate oxygen contents, there are regions of the flake that are highly functionalized as well as regions with no oxygen-functionalities producing a GO flake with a rough surface. At high oxygen contents (C-to-O ratio = 2.0), most carbon atoms are either functionalized themselves or adjacent to at least one functionalized carbon. As a result, the flake has a more consistent thickness and its roughness decreases.

In an aqueous separation experiment, the ability of a solute to pass through the membrane will depend on the connectivity of the water in the nanopores. The connectivity of water was assessed using a recursive algorithm by identifying distinct clusters of water in the membrane models based on a simple distance criterion. A pair of water molecules were defined as being as in the same cluster if the distance between their water oxygen atoms was less than a specified cutoff distance or if they could be connected to each other via intermediary water molecules within that same cutoff. The cutoff (0.36 nm) was defined as the position of the first minimum in the O–O pair distribution function (Figure S5), obtained from a separate simulation of bulk water. For each model, the cluster with the largest number of water molecules was identified and expressed as a fraction, fw, of the total number of water molecules in the membrane model. Thus, fw = 1 indicates that all water molecules in the membrane are in the same cluster and the membrane water is fully connected, and fw ≈ 0 means that all water molecules are isolated or in very small clusters and the membrane has poor connectivity. Figure shows the change in fw with membrane water content and flake oxygen content.

Figure 5

Fraction of water molecules in the largest cluster for membranes with variable water content (black) and flake oxygen content (red). The clusters shown in panels i–iv are example images of the largest water cluster found for the water contents labeled on the main figure.

The increase in fw with water content suggests water connectivity significant improves upon swelling and hydration of the membrane. At high water content (mw > 20%) the water nanopores are fully connected (fw ≈ 1.0). However, at lower water content (mw < 20%) the connectivity of water decreases (fw = 0.64, 0.27, and 0.07 for mw = 15, 10, and 5%, respectively). Therefore, a sharp cutoff in connectivity is apparent at intermediate water content (mw = 20%), above which all of the water in the membrane is connected and below, which it is not. At low water content, the identity of the largest water cluster does not typically change over the course of the 10 ns simulation, suggesting that the water molecules and clusters themselves are relatively immobile. For each model, the largest water cluster can be visually inspected along the profile of fw, with examples shown in Figure i–iv, demonstrating the improving connectivity of the nanopores with water content. Figure shows all of the clusters containing more than 10 water molecules for a membrane model with mw = 10%.

Figure 6

Water clusters, colored separately, for all clusters containing more than 10 water molecules for a model with a water content of 10%. The water molecules not shown are either isolated or present in clusters smaller than 10 molecules.

Figure also shows that the connectivity of water in the membrane is sensitive to flake oxygen content. For mo = 0%, fw = 0.98, indicating that the water in the membrane is almost fully connected. The connectivity decreases to fw = 0.5–0.7 at high oxygen content. Given that the number of water molecules is constant across this range, the improved connectivity at low mo is due to the change in distribution of water molecules throughout the membrane, i.e., the preference of water to reside at the flake edges rather than in the interlayer space (Figure ). This suggests that both the total oxygen content and the distribution of oxygen functionalities (e.g., plane or edges) could be exploited to control the connectivity of water and morphology of nanopores.

The effect of membrane water content and flake oxygen content on the self-diffusion coefficient of water, Dw, was investigated (Figure ). In general, the results show that the diffusion of water in GO membranes is significantly lower (Dw < 1.8 × 10–5 cm2 s–1) than bulk water (Dw = 5.6 × 10–5 cm2 s–1 for the TIP3P model), contradicting some experimental interpretations[2,55] but in agreement with other molecular simulation studies in the literature.[35,40,47,49] The lower value of Dw can be attributed to confinement effects and the hydrophilic nature of GO, which results in strong hydrogen-bonding interactions between nanopore water molecules and the oxygen functionalities of GO flakes.[48,56,57]Dw increases with water content, from 0.1 × 10–5 cm2 s–1 for mw = 5% to 1.8 × 10–5 cm2 s–1 for mw = 40%. The improvement in Dw is slow initially (mw < 15%) but becomes increasingly rapid at higher water contents. This is because at low water content additional water molecules adsorb to available oxygen functionalities of the GO flakes. The number of occupied functionalities increases with water content, and once it reaches a critical value, all of the binding sites in the GO framework are occupied. Additional water molecules must then occupy sites less strongly bound to the surface that are in a more bulk-like water environment. This interpretation is supported by experimental evidence from neutron scattering[48] and broadband dielectric spectroscopy[58] studies of hydrated graphite oxide. These studies showed that the rotational and translational motion of water at low water content is severely restricted, but at higher water contents, motions due to bulk-like water are observed (i.e., associated with translational motion through the membrane). The water content at which the latter motions were first observed in these studies was mw = 15% to 17%, approximately equal to the water content at which Dw begins to rapidly increase in our models (Figure ). It is also notable that this is the water content at which the network of water in the membrane nanopores becomes fully connected (Figure ).

Figure 7

Water self-diffusion coefficient, Dw, for membranes with variable mw (black, circles) and mo (red, squares). For most data points, the error bars are smaller than the symbols.

Following from the explanation that slow water dynamics is due to strong interactions between water and the GO oxygen- unctionalities, it follows to expect a decrease in Dw with flake oxygen content at constant membrane water content. From Figure , it can be seen that this is initially the case, where Dw decreases from 0.77 × 10–5 cm2 s–1 at mo = 0% to 0.15 × 10–5 cm2 s–1 at mo = 15%. Surprisingly, beyond mo = 20%, the diffusion coefficient starts to increases with oxygen content. This can be explained by the ratio of the number of oxygen atoms in water to the number of oxygen atoms in GO (noGO/now) and by the location of the oxygen functionalities. At mo = 0%, noGO/now < 1, water is predominantly found at the flake edges rather than the interlayer space, which does not contain any oxygen functionalities and is therefore more hydrophobic. This leads to the formation of larger clusters with bulk-like water molecules (Figure a). As the oxygen content is increased and basal plane becomes more hydrophilic, more of the water molecules become strongly bound to oxygen-functionalities in the interlayer space, and Dw decreases. At intermediate oxygen contents (mo = 15% to 20%), noGO/now ≈ 1, meaning that each water molecule can be bound to an oxygen-functionality with few unoccupied binding sites (Figure b). At higher oxygen contents (mo > 20%), noGO/now > 1, meaning that every water molecule can be bound to an oxygen-functionality, but there are additional unoccupied binding sites that facilitate the motion of water molecules throughout the membrane nanopores, leading to an increase in Dw. The models prepared provide a basis for the study of the transport properties of GO membranes using simulations, including the quantification of water and solute permeability. This could be investigated using non-equilibrium MD[35] or the kinetic MC approach recently proposed by Apostololou et al.(59)

Finally, the in-plane elastic modulus, E, of each membrane model was calculated (Figure ). Large uncertainties were observed in the calculation of E due to significant variations between models with the same composition. Despite this, at low water content (mw ≤ 15%), the value of E is within the wide range reported experimentally (E = 15–42 GPa),[1] and a statistically significant decrease in E was observed for higher water contents (mw ≥ 20%). The modulus decreases at the same water content that the water network becomes fully connected (Figure ) and is therefore likely to be related to the formation of distinct layers of water (Figure ). A coarser modeling approach, such as that used to study the poromechanics of multilayered clay materials, is required to investigate other mechanical propeties.[60]

Figure 8

In-plane tensile moduli, E, for membranes with variable membrane water content, (black, circles), mw, and flake oxygen content (red, squares), mo.

Conclusions

Understanding the relationship between membrane composition and structure is critical for optimizing the performance of a given membrane separation. For example, understanding how the interlayer distance, nanopore connectivity, and water diffusivity in GO membranes change with relative humidity must be understood to underpin the design of membranes with optimized water flux and salt rejection.[6] The key findings from this work, relevant to the design and optimization of GO membranes for aqueous separations, can be summarized as follows:

GO membranes have a distribution of interlayer distances that is dependent on the membrane water content. For a water content of 15% (dGO = 0.80 nm), the membrane water is predominantly found in monolayers, but dehydrated regions are also observed. For a water content of 40% (dGO = 1.35 nm), dehydrated regions no longer exist and water monolayers, bilayers, and trilayers coexist.

There is a sharp cutoff in nanopore water connectivity at a water content of 20%, below which there are a number of separate water clusters, suggesting that the permeation of solutes in membranes with water contents less than this cutoff may be significantly hindered.

Water self-diffusion coefficients range from 0.1 × 10–5 to 1.8 × 10–5 cm2 s–1 depending on the water content, but diffusion is always slow relative to bulk solution.

All of the properties investigated were found to depend on the oxygen content of the constituent GO flakes, suggesting that separation performance could also be enhanced by optimizing oxygen content, which could be achieved by oxidation or reduction reactions.

Although the direct link between relative humidity and the resulting membrane structure was not made in this study for reasons of computational expediency, this relationship could be established using molecular simulations with a constant imposed water chemical potential. However, because water uptake into the membrane induces significant swelling, the volume of the simulation cell must be allowed to change, so adsorption simulations in a rigid GO framework using the traditional Grand Canonical Monte Carlo (GCMC) approach would be unsuitable. Possible alternative approaches for investigating adsorption-induced structural transitions in flexible adsorbent materials include simulations in the constant pressure Gibbs ensemble,[61,62] osmotic ensemble[63,64] or performing GCMC adsorption simulations coupled with isothermal−isobaric (NPT) MD relaxation cycles.[65]

The intrinsic heterogeneity of the membrane structure highlights the importance of flexible GO flakes in understanding the swelling of GO membranes. Over-simplified models (e.g., in which GO nanopores modeled as rigid, pristine graphitic plates) would not allow the efficient packing of GO flakes (apart from at very low oxygen content, where the flake surface is very smooth), reproduce the gradual swelling with increased water content, or reproduce the coexistence of hydrated and nonhydrated regions between a given pair of flakes at intermediate water content. Therefore, our GO membrane models can be considered a significant improvement upon those previously published, addressing a recently highlighted key research need.[29] Furthermore, the high-throughput and in silico approach adopted and presented in this study is broadly applicable to other layered or porous materials with intrinsic nonstoichiometry and structural heterogeneity and provides a general strategy for their design and optimization.

Methodology

A full description of the procedure and justification for the models of GO flakes is provided in our previous work.[28] In brief, a pristine graphene flake was converted to GO by adding oxygen-containing functionalities until the specified target oxygen content was reached. The size of the flake, 5 nm2, was chosen on the basis on our previous study, which demonstrated that finite size effects were not significant for flakes with an area larger than this.[28] Flake edges were terminated with hydrogen or carboxylic acid functionalities (selected randomly with equal probability) and the flake plane functionalized with hydroxide or epoxide groups (selected randomly with equal probability), justified by experimental observations.[23,30] The random approach to functionalization ensures that GO flake models with the same composition have a different molecular-level distribution of functionalities (Figure S1). Apart from extremely high or low oxygen content, this procedure leads to an amorphous surface, in general, with some small nanocrystalline regions of pristine (aromatic) graphene, which is supported by experimental evidence.[51,66,67] While there is some debate around the size of these regions, experimental evidence shows that the majority of the GO flake surface is oxidized or amorphous. This property is reproduced by our models. The procedure used to functionalize the flakes was the same irrespective of the oxygen content, and the specific types and distributions of oxygen functionalities were not modified to reflect the distinct surface chemistry of RGO for low oxygen content.[24,26]

Each GO membrane model was constructed according to the following procedure. A total of six GO flakes were randomly distributed, along with water, in a simulation cell with side lengths of L = 5 nm, L = 5 nm, and L = 50 nm, with a structure-less wall placed at L = 0. Interactions between the membrane atoms and the wall were modeled using Lennard-Jones 12–6 potential (σ = 0.355 nm, ε = 0.293 kJ mol–1). The centers of mass of the GO flakes and water molecules were slowly pulled toward the base of the simulation cell during a 2.5 ns NVT simulation, at a rate of 10–3 nm ps–1, using a harmonic potential with a force constant of 103 kJ mol–1 nm–2, whereupon the GO layers deposited on the wall. Following this simulation, the wall was removed and a 10 ns anisotropic NPT simulation was performed to remove remaining vacuum gap, resulting in a periodic membrane structure in 3D. A further 10.1 ns isotropic NPT simulation was performed, from which the time-averaged interlayer distance distribution, water connectivity and water diffusion coefficient were determined from the final 10 ns, using atomic coordinates collected every 200, 100, and 20 ps, respectively. Diffusion coefficients were obtained using the Einstein relation, from the gradient of a plot of mean-squared displacement against time at short time scales (between 5 and 25 ps). Pore size distributions were determined from the final configuration after removing membrane water molecules using the Poreblazer software.[68] In-plane tensile moduli, E, were calculated from a straight line fit to the stress–strain relationship in the elastic region (<0.5% strain), obtained from additional 10 ns simulations in which the simulation cell is gradually extended in the x direction at increments of 1.0 × 10–5 nm every 1 ps. In these simulations, the cell lengths in the y and z directions are allowed to relax in response to the deformation.

All MD simulations were performed using GROMACS, version 2016.3.[69] GO and water were modeled using the CHARMM[70] and TIP3P[71] force fields, respectively. van der Waals interactions were modeled using a LJ 12–6 potential that was smoothly cut off between 1.0 and 1.2 nm and parameters for atoms that were not alike were obtained using the Lorentz–Berthelot combining rules. Equations of motion were integrated using the leapfrog method[72] using a time step of 1 fs, and TIP3P water molecules were constrained using the SETTLE algorithm.[73] Electrostatic interactions were evaluated to within a tolerance of 10–5 kJ mol–1 using the PME approach[74,75] with a real-space cut-off of 1.2 nm. A target temperature of 298 K was maintained using the Nosé–Hoover thermostat[76,77] with a coupling time constant of 1 fs. In the constant pressure simulations, a target pressure of 1 bar was maintained using the Parrinello–Rahman barostat,[78] with a time constant of 1 fs and a compressibility of 4.5 × 10–5 bar–1.