Molecular dynamics simulation of water transport through graphene-based nanopores:Flow behavior and structure characteristics☆

Xueping Yang ,Xiaoning Yang ,*,Shuyan Liu

1 State Key Laboratory of Material-Orientated Chemical Engineering,College of Chemistry and Chemical Engineering,Nanjing Tech University,Nanjing 210009,China

2 College of Mechanics and Power Engineering,Nanjing Tech University,Nanjing 210009,China

Keywords:Graphene Nanopores Flow behavior Membrane Molecular simulation

ABSTRACT The flow behavior ofpressure-driven water in filtration through graphene-based slitnanopores has been studied by molecular simulation.The simulated flow rate is close to the experimental values,which demonstrates the reasonability of simulation results.Water molecules can spontaneously in filtrate into the nanopores,but an external driving force is generally required to pass through the whole pores.The exit of nanopore has a large obstruction on the water effusion.The flow velocity within the graphene nanochannels does not display monotonous dependence upon the pore width,indicating that the flow is related to the microscopic structures of water con fined in the nanopores.Extensive structures of con fined water are characterized in order to understand the flow behavior.This simulation improves the understanding of graphene-based nano fluidics,which helps in developing a new type of membrane separation technique.

1.Introduction

Graphene-based materials,such as graphene and graphene oxide,have been considered as promising membrane materials[1-4].Graphene layers usually self-assemble into laminate paper-like structures with interlayer distance on nanometer scales[5-7].This speci fic laminate structure allows watermolecules to permeate through interconnected nanochannels between graphene nanosheets.Recently,Nair et al.[3]have reported submicrometer-thick graphene-oxide laminate membrane with a typical～1 nm pore width.This kind of membrane can impede the permeation of various species,including helium,but the membrane allows unhampered permeation for water.As further extension of the pioneer work,Joshi et al.[8]evaluated the filtration and separation performance of laminate membrane for extensive solutes,including ions and organics.In addition,ultrathin graphene nano filtration membrane with two-dimensionalnanocapillary has been fabricated by packing reduced graphene oxide structures[4].The graphene membrane performs excellently for retention of organic dyes with high water flux.

Although carbon nanotubes(CNTs)have been reported as potential nano filtration membranes,the CNT-based membranes suffer from the dif ficulty in larger-scale fabrication and preparation.On the contrary,the laminate two-dimensional graphene membranes,with superior flexibility and chemical stability,can be facilely prepared by filtrationassisted assembly method.According to a previous model[3],the interconnected nanochannels in the layered graphene membrane have two regions:functional and pristine.The former acts as spacers to keep adjacent graphene sheets apart,whereas the pristine region provides a capillary network that allows high water flux.It is also suggested that the nanoscale spacing between graphene nanosheets can be adjusted through inserting spacers[9].Therefore,a broad spectrum ofgraphene membrane could be prepared,which shows extensive potential applications.

In order to develop this new-typed graphene-based laminate membrane,it is necessary to understand the flow behavior and con fined structure of fluid molecules in the two-dimensional nanochannels.Molecular simulation has been extensively applied to study the transport behavior and con fined structures in the nanoscale CNTs[10-12].However,to our best knowledge,no molecular simulation has been reported for the pressure-driven in filtration transport of water molecules from bulk phase going through the slit nanochannels formed by graphene sheets.In particular,graphene pores with subnanometer size show signi ficant molecular sieve performance,which can be used to develop new-typed membrane separation materials[3,9].For this kind of pore size,the interfacial effect is expected to overwhelmingly control the flow behavior.It is highly required to explore the unique flow behavior through the subnanometer graphene slit pores.In this work,we use molecular dynamics(MD)simulation to study the in filtration behavior of water driven by external pressure through the graphene nanocapillary with typical pore widths(≤1 nm).The con fined water structures are also investigated in order to understand the flow mechanisms.

2.Simulation Methods

In this work,water molecules are simulated using the SPC/E model[13].Carbon atoms of graphene are held stationary and modeled as Lennard-Jones(LJ)spheres employing the parameters proposed by Chang and Steele[14].The van der Waals interactions between different particles are calculated with the LJpotentialusing the Lorentz-Berthelot mixing rule.This water-carbon interaction has been successfully applied to represent the structural and dynamic properties of water molecules con fined in carbon pores[15].

As shown in Fig.1,the graphene-based slit nanopores consist of a pair of three-layer graphene sheets,which adopt ABA stacking con figurations[16]with the natural interspacing space of 0.335 nm[17].Each graphene sheet has the size of 2.333×5.105 nm2.The slit pore is incorporated into the simulation cell with various separation distances(H=0.7,0.8,0.9,and 1 nm).The simulation box lengths in x,y,and z directions are set to be 2.48,2.40,and 17.21 nm,respectively.Water molecules are randomly placed in the right reservoir with a density of 1.0 g·m-3and the leftreservoiris keptempty.One additionalgraphene sheet is as moveable wall for creating the driven force toward the water molecules in the right reservoir.The simulated pressure of the right water reservoir is obtained around 0.1 MPa.

We perform equilibrium and non-equilibrium MD simulations in the canonical ensemble(NVT)using the LAMMPS MD package[18]with a time step of 1 fs.The Nose-Hoover thermostat[19,20]is used to keep the temperature of 298.15 K.The particle-particle particle-mesh method is used for the Coulombic interaction.The cutoff for LJ potential is set to be 1.2 nm.Periodic boundary conditions are used only in the x and y directions(the coordinate system in Fig.1).The con fined structures are mainly evaluated through the equilibrium MD simulation.In the non-equilibrium MD simulation,following the general treatment in Ref[11,12],an externalforce along the z direction isapplied to the moveable wall in order to push water molecules to enter the nanopores.The applied pressure(P)is calculated from the moving acceleration rate(a)of moving graphene sheet,P=Nma/A,where N is the number of atom in the con fined pore,m is the mass,and A is the area ofmoving graphene.In this work,we use relatively high pressure,which is common in the non-equilibrium MD simulation[12]to reduce the thermal noise and speed up the simulation process.Meanwhile,in order to simulate the water penetration in the pervaporation membrane process,the water molecules after passing the slit pore will be removed from the system.

3.Results and Discussion

Fig.2(a)shows thatthe volumetric flow rate(Q,nm3·ns-1)ofwater through the slit nanopores increases with the applied pressure.Fig.2(b)shows an approximate linear dependence of the number of permeating molecules on the time,fromwhich the flow rate is obtained fromslopes ofthese curves.Moreover,effective flow velocity is estimated to be in the range of(1.5-27.8)×104nm·s-1in the pore of 1 nm,which is higher than the experimentally observable flow velocity of water within the CNT pore(～0.43 × 103nm·s-1).It should be noted that under low-pressure condition,the linear dependence behavior is somewhat suppressed due to thermal noise effect.As expected,larger slit nanopore leads to higher volumetric flow rate.It is noteworthy that the volumetric flow rate does not exhibit a linear relationship with pressure.In the graphene nanopores with the pore width of 0.7-1 nm,the simulated mass flow rate varies from 1.9×10-11to 6.4 × 10-10mg·s-1in the pressure ranges investigated,which are consistent with the experimentally measured permeability of water through graphene oxide membranes with the width of 0.7-1.1 nm[3].The result suggests that this simulation system can describe the basic flow characteristics of fluid through graphene-based slit nanopores.

We further evaluate the flow enhancement factor(ε=Q/QHP)in our graphene nanopores,in which the Hagen-Poiseuille equation[21,22]is used to describe the classical flow rate(QHP).The inset of Fig.2(a)displays the enhancement factors,which are lower than those reported for water flow through CNTs[21,23].The enhanced flow behavior is generally explained by small friction between water molecules and carbon surfaces,corresponding to the so-called boundary slip behavior[24-26].However,different carbon-based nanopores often exhibit variant flow features.The simple boundary slip concept cannot give a full explanation to the microscopic characteristics of fluid flow through the graphene-based slit nanopores.Especially,in our system,there exist capillary suction and discharge between bulk phase and micropore phase,which might affect the flow velocity.

The simulation results indicate that water molecules could spontaneously permeate into the microscopic nanopores,because the capillary pressure generated at the pore entrance would drive water into the nanopores.Without the applied pressure as the driven force,water could not permeate the nanopores with the pore width below 10 ?.Moreover,Fig.2(b)shows that at the initial stage water molecules could not enter the nanopores under lower pressure.The reason is that certain external force is required to overcome the permeation resistance from micropore surface interaction[3,27,28].In our simulation system,the interaction between graphene surface and water molecules exhibits somewhat hydrophilic characteristics[29],as a result,the capillary permeation pressure at the exit of nanopore would block or suppress fluid permeation.

Fig.3(a)shows the variation in the number of water molecules entering the nanopore(H=0.8 nm)versus time.Higher pressure results in larger driving force,increasing the transportation rate.Similar phenomena are also observed in other size nanopores.The permeation and flow behavior of water molecules in the graphene-based slit nanopores involve three steps:(1)permeate into nanopores;(2)move along nanopores;and(3) flow out of nanopores(exit effect).Fig.3(b)shows the axial(z-direction)trajectories of several typical water molecules through the 0.7 nm pore under the pressure of 108 MPa.In the graphene slit nanopores,the axial movement of water molecules demonstrates a typical pulse-like behavior to some extent,which is in agreement with the behavior of water molecules in singlewalled carbon nanotubes[30].However,as compared to the moving behavior along the direction parallel to the surface,the movement in the lateral direction does not show signi ficant fluctuation,demonstrating the restriction effect of solid surface.Moreover,Fig.3(b)shows that the dwelling time ofwatermolecules near the nanopore end is obviously long,indicating that the permeation of water molecules is largely hindered by the exit effect.Generally,this exit effect can be quantitatively characterized by the free energy pro file for a molecule passing through the exit[31],which re flects the degree of diffusivity for molecule passage.However,in this work,this computation is not conducted.

Fig.1.Snapshot of system schematic in the non-equilibrium MD simulation.

Fig.2.(a)Volumetric flow rate(Q)versus pressure,with the insetshowing the variation in flow enhancementfactor(ε)with pressure;(b)the numberofwater molecules permeating the graphene slit nanopores as a function of simulation time.

Fig.3(c)shows the average flow velocity ofwatermolecules through the nanopores.This velocity,ranging from 5 to 20 m·s-1,is consistent with previous results[24,32].The diverse velocities for different pore widths re flect the flow resistance effect[21].Clearly,the flow velocities increase with the pressure and do notshow monotonous increase ordecrease trend with the pore width.At H=0.8 nm the flow velocity is the largest,which can be correlated with the con fined structures(see the nextsection).The diffusion mobility ofwater molecules in the graphene slit pores is also characterized by the parallel mean square displacements(MSD||)along the flow direction,which is computed according to

where the quantity in the broken bracket averaged over many different initial time t0and x and z represents the molecular positions in the two spatial directions.

As shown in Fig.3(d),the MSD||is substantially less than that of bulk water,indicating the pore con finement effect.The inset of Fig.3(d)gives the diffusion coef ficients ofwater in differentsizes of nanopores.The dependence of diffusion coef ficients on nanopore width is in accord with the relationship between the flow velocity and pore width,thatis,higher diffusion coef ficient for H=0.8 nm corresponds to a larger flow velocity(see Fig.3(c)).This observation means that the con fined water structure plays a critical role in the transportation of water through the graphenebased slit nanopores.This con fined structure will be evaluated in detail in the following section.In addition,as shown in Fig.3(c),in the pore of H=0.7 nm there is a jump in the flow velocity as the pressure increases,which is probably ascribed to the structural change in the con fined water molecules.

Fig.3.(a)The number of water molecules inside slit nanopore(H=0.8 nm)as a function of time;(b)moving trajectories of several water molecules through the slit pore(H=0.7 nm)with the pressure of 108 MPa;(c)the flow velocity of water along the z-axis direction as a function of pressure;(d)the computed MSD||in different pore con finements,with the inset showing the diffusion coef ficient of water molecules.

Fig.4.(a)The water density pro file along the direction perpendicular to graphene surfaces under the pressure of 108 MPa;(b)the effective water density(ρe)as a function of nanopore width;(c)the lateral radial distribution function g o-o(r)for water layers adjacent to the pore surfaces;(d)the snapshots of monolayer water structures inside nanopore with H=0.7 nm.

In order to further understand the water penetration behavior,we simulate and characterize the con fined structures for water flow through the graphene-based nanopores.Fig.4(a)shows typical density pro files under the pressure of 108 MPa.Similar result can be observed under other pressures.The effective water density is given in Fig.4(b).In general,the density pro files are symmetric with respect to the pore center.Single layer of water molecules forms in the slit width for H=0.7 nm and two layers at H=0.9 and 1 nm.A transition layer between single and double layers appears for H=0.8 nm(see Fig.4(a)).The results agree with previous simulations for water con fined inside graphene slitpores[7,17,33].However,the density pro file isusually associated with its con fined density.In the previous work[34],it has been reported that when the con fined water density decreases from 1.48 to 1.11 g·cm-3,the density pro file merges from two water layers into single layer.In this work,the con fined water density at H=0.7 nm is approximately 1.0 g·cm-3,so the formation of single layer density pro file is reasonable.

There is interfacial density depletion region near the graphene surfaces in the density pro files(see Fig.4(a)).This might correspond to the interface boundary slip behavior,leading to the rapid water flow in the graphene pores.From Fig.4(b),the con fined water density follows the order ρe(H=7)＞ ρe(H=9)＞ ρe(H=10)＞ ρe(H=8).As expected,higher con fined density leads to slower diffusion of water molecules.This has been con firmed by the neutron scattering experiment[35].A careful inspection of Fig.4(b)indicates that the density slightly reduces with the increase of pressure.This observation is consistent with the minor pressure dependence of con fined water density in CNTs[36].

Fig.4(c)presents the radial distribution functions(RDFs,go-o(r))of con fined water molecules within the surface layer along the direction parallel the surface at two pressures.We observe the first peak of go-o(r)at～0.27 nm,representing the first coordination shell between water molecules.This is consistent with the structural behavior of bulk water.The RDFs of water con fined in the graphene pores agree well with previous results[7,34].As the pore width decreases,the second peak position in the con fined RDFs shifts toward large separation.In particular,for H=0.7 nm,the con fined RDF(red curve)has the characteristic of a solid-like oscillation,exhibiting a well-de fined long-range order structure.This RDF behavior is related with the 4-membered ring structures of con fined water molecules in the 0.7 nm pore,as given in the snapshots(see Fig.4(d)).According to our simulation result,the graphene nanopores can induce ordered con fined water structures,which is highly connected with the anomalously fast movement of water molecules in the nanopores.Our result provides a further support on the previous speculation that an ordered structure forms as water flows in graphene laminate membranes[3].It is interesting to note that for H=0.7 nm,as the pressure increases,the 4-membered water-ring structure has certain distortion.The con figuration distortion can be further con firmed by the reduced ratio between the peak and valley heights in the RDF curve of the 0.7 nm pore under higher pressure(red line in Fig.4(c)).This may be the reason why the flow velocity in the pore of H=0.7 nm has an obvious increase with the pressure.

Fig.5 gives the orientationalbehaviorofcon fined water molecules in the graphene-based slit nanopores.Two types of orientational angles(de fined in Fig.5)are analyzed.The result shows that water molecules mainly orient themselves parallel to the graphene surface,which is consistent with the previous result of the water orientation near graphene surface[15].This parallel orientational behavior is in favor of the formation of hydrogen bonds(H-bonds)among the con fined water molecules.Fig.6(a)shows the hydrogen-bond pro file of the water molecules in the graphene slit pores.A geometric criterion is used to de fine the hydrogen bond[37,38].The H-bond pro files are symmetric with respect to the pore center and the H-bond number is near 3.5 in the center region.This H-bond behavior is in line with the previous simulation[39]and the vibrational spectroscopy[40]of interfacialwater H-bonds.The average number of H-bonds is from 3.1 to 3.5,depending on the pore width.For H=0.7 nm,the number ofH-bonds is the largestand itis the leastfor H=0.8 nm.Itis already known thatenhanced H-bond structure could suppress the mobility of water molecules inside slit nanopores[27,34].In addition,in the 0.8 nm pore,the con fined density of water is the least,along with weak interfacial structure layer,as shown in Fig.4(a,b).This willproduce low flow resistance.On the basis of these analyses,the water movement should be the fastest inside the nanopore of 0.8 nm,demonstrating larger flow velocity through the nanochannel.This is consistent with the previous simulation result.

Fig.5.Orientationaldistributions of water molecules in slitnanopores.(a)OHbond orientations;(b)dipole moment orientations.

The H-bond stability for con fined water molecules can be characterized by the intermittent time correlation function,CHB(t),which is de fined as[41,42]:CHB(t)= 〈h(0)h(t)〉/〈h(0)h(0)〉,where h(t)is unity when a particular H-bond pair inside nanopore exists at time t,and zero otherwise.The calculated CHB(t)function for the H2O-H2O H-bonds inside the slit nanopore is shown in Fig.6(b).For comparison,the bulk CHB(t)function is also given.The CHB(t)curves of water molecules inside slit nanopores decay slower than the bulk CHB(t),meaning that the H-bond strength of con fined water is enhanced.We note that the decay rate of the correlation function is the slowest for H=0.7 nm.This can be explained in terms ofthe enhanced hydrogen bond interactions caused by the 4-membered water-ring structure.Meanwhile,as can be seen from Fig.6(b),the decay rate of CHB(t)would increase as the slit width increases,so the stability of H-bonds will become weak and the flow rate of water in the slit pore will increase.

4.Conclusions

In this work,MD simulations are performed to study the in filtration behavior of pressure-driven water flow through two-dimensional graphene nanochannels with the widths from 0.7 to 1 nm.Simulated flow rates are close to the relevant experimental data.This demonstrates that the simulated system could reasonably re flect the flow characteristics of the graphene-based slit nanopores.The flow enhancementof the water through microscopic nanopores is smaller than those reported for CNTs.There exists pulse-like transmission of water penetration through nanopores and the in filtration of water molecules through slit nanopores is slowed down by the pore exit.The flow velocity of water molecules through the nanopores does not display monotonous dependence on the pore width.The results illustrate that the flow velocity is related to microscopic structure of con fined water.The extensive structure analysis shows that con fined water molecules form ordered layer structure,which is expected to be favorable to the fast flow of water through the graphene-based slit nanopores.The flow velocities for various pore widths can be interpreted in terms of the con fined structures and hydrogen bond network.

Fig.6.(a)Average number pro file ofhydrogen bond for water along the direction perpendicularto graphene surfaces;(b)the intermittenttime correlation function C HB(t)for H2O-H2OH-bonds.