- Email: [email protected]
Transport of salty water through graphene bilayer in an electric field: A molecular dynamics study
Computational Materials Science 131 (2017) 100–107
Contents lists available at ScienceDirect
Computational Materials Science journal homepage: www.elsevier.com/locate/commatsci
Transport of salty water through graphene bilayer in an electric field: A molecular dynamics study Hui Zhang a, Bo Liu a, Mao-See Wu b, Kun Zhou a,b,⇑, Adrian Wing-Keung Law a,c,⇑ a Environmental Process Modelling Centre, Nanyang Environment & Water Research Institute (NEWRI), Nanyang Technological University, 1 Cleantech Loop, CleanTech One, 637141, Singapore b School of Mechanical and Aerospace Engineering, Nanyang Technological University, 50 Nanyang Avenue, Singapore 639798, Singapore c School of Civil and Environmental Engineering, Nanyang Technological University, 50 Nanyang Avenue, Singapore 639798, Singapore
a r t i c l e
i n f o
Article history: Received 1 November 2016 Received in revised form 24 January 2017 Accepted 28 January 2017
Keywords: Molecular dynamics Graphene bilayer Electric field Salt rejection
a b s t r a c t Molecular dynamics simulations are carried out to study the pressure-driven transport of salty (NaCl) water through the nanochannels formed by a graphene (GE) bilayer with and without a vertical electric field. The simulation results show that the channel thickness influences not only the velocity and numbers of water molecules and ions inside the channel, but also their responses to varied pressure and electric field intensity. The salt rejection rate is characterized by considering both the numbers of ions and water molecules inside the channel and their velocities. When the electric field is not imposed, the salt rejection rate varies slightly in the high-pressure range (50–300 MPa), and is mainly affected by the channel thickness. When an electric field with high intensity (1–5 V/Å) is imposed vertically under high driving pressure (300 MPa), ions can be trapped inside the nanochannel. Simulations under lower electric field (0.01–0.1 V/Å and 0.05–0.2 V/Å) and pressure (5 MPa and 10 MPa) are also conducted. Under these circumstances, GE nanochannels can reject ions from entering the channel, indicating improved salt rejection rate and reduced ion channel accumulation. Ó 2017 Elsevier B.V. All rights reserved.
1. Introduction Graphene (GE) has recently attracted substantial attention in water desalination techniques [1–8]. Reverse osmosis (RO) membranes derived from GE have been reported to show advantages over conventional RO membranes because of ultimate thinness, high mechanical strength, and more importantly, high water permeability [2–5,9–11]. GE derived membranes can be fabricated by stacking functionalized GE nanosheets, where nanochannel and nanopores are formed. Interlayer distance of graphene oxide (GO) membrane can be less than 0.7 nm if made by vacuum filtration in the dry state, enabling water molecules to permeate through while blocking hydrated ions and other target molecules [5,3,12]. Confined in the nanostructures, water still shows high permeability through GO membranes. This is attributed to the capillary network formed by the pristine GE domains, where water can flow almost frictionlessly [5,11].
⇑ Corresponding authors at: Environmental Process Modelling Centre, Nanyang Environment & Water Research Institute (NEWRI), Nanyang Technological University, 1 Cleantech Loop, CleanTech One, 637141, Singapore. E-mail addresses: [email protected] (K. Zhou), [email protected] (A.W.-K. Law). http://dx.doi.org/10.1016/j.commatsci.2017.01.039 0927-0256/Ó 2017 Elsevier B.V. All rights reserved.
GE-based membranes also have superior energy efficiency compared with other RO technology [13]. Meanwhile, functionalized GO can enhance the hydrophilicity of the membrane and provide negative charges, and thus possessing anti-fouling properties [14]. However, the scale-up ability remains a critical issue for the mass production of GE-based membranes. Many studies have reported on various possible fabrication approaches, such as chemical vapour deposition and subsequent transfer to substrates [15], restacking GO flakes by filtration on a back support [2], etc. Most of these methods suffer from scalability issues. Recently, Akbari et al. developed a method to produce large-area, highly ordered, and multi-layered GE-based films by making use of the shear alignment of discotic nematic GO fluid [16]. This scalable method brings hope to the mass production of GE-based membranes for industrial applications. Molecular dynamics (MD) simulation has been used to investigate the transport of salty water through GE nanopores and nanochannels [4,11,17–24]. GE single layer with nanopores bonded with chemical functional groups at the pore ring can limit the salt passage depending on the pore size, functional groups and applied pressure or electric field [18–22]. Cohen-Tanugi et al. also carried out MD simulations of salty water passing through multilayer nanoporous GE with different interlayer width as well
101
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
as pore offset, and demonstrated that salt separation performance can be manipulated by different multilayer design [17]. Simulation results also showed that the water molecules display layered structure inside GE nanochannels which affects the water flow rate significantly [4,8]. GE nanochannels with the thickness of allowing one monolayer of water molecules passing through do not allow ions inside [4]. However, a systematic study of ion transport through the nanochannel formed by GE bilayers, especially the interaction of ions with water molecules and GE, is still needed to provide further guidance on the design and application of GEbased RO membranes. High surface area and excellent conductivity also make GE feasible for the fabrication of electrodes in capacitive deionization (CDI) [1,8,26,27]. The experimental results from these studies showed that the electrosorption capacity of GE-derived electrodes can be enhanced by mixing GE with other chemicals or creating specific structures. In CDI, salty water flows between the electrodes [26,27]. MD simulations can provide insights into the effect of nanostructure on CDI performance which is still not fully determined or optimized [28]. The present work focuses on the behavior of ions and water molecules in flowing solutions rather than static solutions, while studies so far often use the latter one for MD simulations due to their relative ease to be set up [25,29]. Nasrabadi and Foroutan carried out MD simulations to demonstrate that carbon nanotube (CNT) electrodes are suitable for ion separation by investigating the atomic interaction among ions, water molecules and CNT [25]. However, the feasibility study of GE derived electrodes is mostly based on material characterization and experiments so far. The present study is targeted at a systematic investigation of the transport of salty water through the nanochannel formed by a GE bilayer with different channel thickness, applied pressure as well as external electric field using MD simulations. Interactions among the ions, water molecules and GE are also investigated. The simulation results provide useful guidelines for the design and application of GE-derived material for RO or CDI desalination.
2. Simulation details The studied system consists of a nanochannel formed by the bilayer GE and two water reservoirs at each end of the channel constrained by GE walls, as shown in Fig. 1. The left-hand feeding
reservoir is initially filled with salty water containing 16 Na+ and 16 Cl dissolved with 922 water molecules, and the right-hand receiving reservoir with pure water containing 954 water molecules. The salinity in the feeding reservoir is 4.86 g/L, slightly higher than sea water salinity (3.5 g/L). Both reservoirs have the same dimension of 40 20 40 Å3. The GE bilayer is parallel to the xy plane with a constant length L 32 Å, while its thickness d varies from 7 to 20 Å. The GE constrain walls are set normal to the bilayer plane to restrict water and ion movement. Simulations are carried out using the large-scale atomic/molecular massively parallel simulator (LAMMPS) package [30]. Boundary conditions are set as periodic in all three dimensions to maintain the continuous feeding. The Transferable Intermolecular Potential 3P (TIT3P) [31] is applied to model the water molecules, and partial charges of O and H are 0.834e and 0.417e. Other atomic interactions are modelled using CHARMM force field [32], including Coulombics and van der Waals (vdW) interactions. The long-range Coulombics are computed with the particle-particleparticle-mesh (PPPM) method with a cutoff of 12 Å, and vdW interaction is modelled using Lennard-Jones (LJ) potential with a cutoff at 10 Å. LJ parameters and charges of each atom species are summarized in Table 1 according to previous research studies [23,32,33]. The LJ pair parameters between species i and j are calculated using the Lorentz-Berthelot combining rules except for C-O and C-H interactions. The energy well e and zero-across distance r between C and O are set as 0.0937 kcal/mol and 3.19 Å, while e and r between C and H are all set as 0 [23]. The NVT ensemble with a Nosé-hoover thermostat is used to keep the temperature at 300 K. The GE nanochannel and constrain walls are treated as rigid body with carbon atoms fixed. The SHAKE algorithm is used to constrain the bond and angle of water molecules to maintain the rigid molecule structure. After all the water molecules are relaxed for 50 ps with a timestep of 0.5 fs to reach equilibrium state, a pressure difference DP is applied by adding a force f on all the O atoms in the left feeding reservoir along the positive X-direction. The force f is calculated as [23]
f ¼ DPLy Lz =n
ð1Þ
where n is the number of water molecules in the left water reservoir, and LyLz is the cross-section area shown in Fig. 1. Our simulations focus on the pressure difference DP ranging from 50 to 300 MPa to obtain data for statistical calculations. The external electric field is added in the channel region and vertical to the xy plane. Because ions and water molecules move very fast due to the high pressure, the applied electric field E also has a large intensity (1–5 V/Å) to obtain distinguishable results. The intensity range here is with the same magnitude of electric field that can break bonds and induce ionization of water molecules [34], but the rigid structure of water molecules is kept using the SHAKE algorithm. Cases with lower pressures (5 MPa and 10 MPa) and lower electric field intensities (0.01–0.1 V/Å and 0.05–0.20 V/Å) are also tested, which are discussed in Section 3.3. The simulated system can reach the steady state within around 100 ps with the imposed force f and electric field E. The simulation is then carried on for another 1.5 ns (2.5 ns when DP = 5 and 10 MPa) to collect data. Table 1 LJ parameters and charges used in the simulation.
Fig. 1. Front and top view of the simulated system consisting of two reservoirs and a nanochannel formed by the GE bilayer. At the initial condition, the left-hand reservoir is filled with salt (Na+ and Cl) water, and the right-hand reservoir is filled with pure water. The two reservoirs are separated by GE nanochannel and constrain wall.
Species
e (kcal/mol)
r (Å)
q (e)
H O C Na+ Cl
0.0461 0.1521 0.0700 0.0874 0.0356
0.40 3.15 3.98 2.74 5.03
0.417 0.834 0 1 1
102
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
3. Results and discussion 3.1. Influence of pressure and channel thickness of In each simulation, the number N and average velocity v water molecules and ions are obtained after the system reaches remain constant with slight fluctuthe steady state where N and v of water molecules and ions are compared ations. The N and v under the applied pressure difference DP varying from 50 to 300 MPa at the channel thickness d of 9 Å, 13 Å and 20 Å, as shown in Fig. 2. The three thicknesses chosen here correspond with twolayer, three-layer and bulk water distribution structures inside the nanochannel, reported by previous studies [4,23,35]. The one-layer water flow at smaller thickness is not chosen because ions can hardly enter such narrow channel [4,17]. An almost linear increase of velocity in response to pressure DP can be observed for both water molecules and ions in Fig. 2 when the electric field is not imposed. At the same channel thickness, the velocities of ions have similar magnitude to but are slightly smaller than the velocities of water molecules. The slight difference of velocities between the water molecules and ions becomes more apparent when the channel thickness is smaller. At a larger channel thickness, the velocities of both water molecules and ions increase more rapidly with the increasing pressure, suggesting that a weaker constraint by the nanochannel enables more significant influence of the driving pressure. The velocity of water molecules inside the CNTs with diameters of 11–16.6 Å at DP = 75 MPa varies from 2 to 6.5 m/s [36]. In the present simulations, the velocities of water inside the GE nanochannel of d = 13 Å are 10.6 m/s and 16.6 m/s at DP = 50 MPa and DP = 100 MPa respectively, higher than that in the CNTs. It can be observed from Fig. 2 that the number of water molecules inside the channel varies little in the pressure range of 50– 300 MPa. As for the number of ions, at d = 20 Å, no obvious change
occurs. At d = 13 Å, the lower pressure of 50 MPa starts to show its effect on reducing N (ION); at d = 9 Å, a more obvious drop of N (ION) can be observed. This suggests that the number of ions inside the nanochannel is more sensitive to pressure changes than water molecules, especially at smaller channel thickness. Dissolved in water, ions are surrounded with hydration shell and they tend to lose part of their hydration shell when confined in other nanostructures [37,38]. In an aqueous NaCl solution of similar salinity with that in this study, the hydration radii of Cl and Na+ are 3.23 Å and 2.35 Å respectively [39,40]. The equilibrium vdW distance between water molecules on the hydration shell and the GE layers is 3.19 Å (equal to rC-O), as illustrated in Section 2. Therefore, in order for ions to get into the narrow nanochannel at d = 9 Å, their hydration shell has to be partly peeled off, which demands more energy. The dependence of transport properties of water and ions inside the channel on its thickness d is described in Fig. 3 with a constant pressure difference and without electric field. Our previous work has reported the dependence of N (H2O) on the channel thickness for d = 7–13 Å in detail, including the sharp increase between d = 8 Å and 9 Å as well as the less sharp increase between d = 11 Å and 12 Å [23]. In this study, d = 20 Å is also investigated to indicate the dependence of d over a wider range. For d = 7 Å, ions cannot pass through. This can be explained by the radii of Na+ and Cl (1.02 Å and 1.81 Å) [41] and their equilibrium vdW distance with the GE layer (rNa-C = 3.36 Å and rCl-C = 4.50 Å). The increase in gradient of N (ION) is similar to that of N (H2O), except there is a sharper increase from d = 8 Å to d = 9 Å for N (ION), where water molecules change from two-layer to threelayer. From the velocity profile, it can be observed that although v (H2O) and v (ION) increase with the similar trend as d gets larger, v (H2O) is always larger than v (ION). The velocity difference indicates that the constraining effect of the nanochannel on ions is stronger than that on water molecules. The difference between
Fig. 2. Influence of pressure on (a) velocity of water molecules, (b) ratio of velocity of ions (average of Na+ and Cl) vs velocity of water molecules, (c) number of water molecules, and (d) number of ions (sum of Na+ and Cl) inside the GE nanochannel of thickness d = 9, 13 and 20 Å when electric field is not applied.
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
103
Fig. 3. Influence of channel thickness on (a) number of water molecules inside the channel, (b) number of ions inside the channel, and (c) velocity of water molecules and ions (average of Na+ and Cl) inside the channel at constant pressure difference DP = 300 MPa.
v (H2O) and v (ION) reduces as d gets larger, because a nanochannel with larger d has a weaker constraining effect.
channel. Therefore, in the present study, the salt rejection SR is defined as below,
3.2. Salt rejection
SR ¼ ð1 KC=C 0 Þ 100%
Salt rejection rate (SR) is an important parameter to assess the performance of RO membrane. In engineering, SR is defined as one minus the ratio of permeate salinity versus feed salinity. In MD simulations of either GE or CNT, the definition of SR generally follows the engineering approach. Previous research calculates SR by comparing the number of permeated water molecules and ions to the initial values in the feed reservoir, or only by counting the number of permeated ions, regardless of whether that boundary condition is periodic or not [17,42–44]. Periodic boundary condition enables the water molecules and ions to move from one reservoir to another without passing through the membrane. Therefore, it brings limitations to the general definition of SR based on counting numbers, which is also highlighted by Thomas and Corry [44]. However, under non-periodic boundary condition, another problem arises regarding when to count the numbers. This is because ion concentrations under non-periodic boundary conditions change from time to time rather than reaching a steady state. Thus, some researches use the time point instead when certain amount of water molecules or ions enter into the permeate reservoir [17,43]. As discussed in Section 3.1, the GE nanochannel can influence the ions passing through by reducing their velocity inside the channel as well as the opportunity for them to enter the channel. These two aspects can be expressed by the comparison of velocity and number between the ions and the water molecules inside the
where K ¼ v ðIONÞ=v ðH2 OÞ is the ratio of average velocities of ions and water molecules inside the channel; C = N (ION)/N (H2O), the ion concentration inside the channel, expressed as the ratio of the average numbers of ions and water molecules inside the channel; and C 0 ¼ N0 ðIONÞ=N 0 ðH2OÞ, the initial ion concentration expressed as the ratio of the initial numbers of ions and water molecules in the left feeding reservoir. The definition of SR here can overcome the difficulties of counting numbers under the periodic boundary condition and help decide on the time point under non-periodic boundary condition. Each simulation first runs for a rather long time so that the system reaches the steady state, simulation data are then collected to calculate the average values of each parameter in Eq. (2) to make SR more statistically reasonable. In fact, given any time period, Eq. (2) expresses the same meaning as that of the general definition by integration over time. In this manner, we can ensure that Eq. (2) calculates a steady-state SR. Indeed, we are not able to obtain large salinity difference at each side of the channel, but under high pressure range, its influence on diffusion could be neglected. Fig. 4 shows the salt rejection under the influence of pressure difference DP and channel thickness d. Within the investigated range of DP, SR at d = 13 and 20 Å under different pressure does not vary significantly, implying that pressure has not started to show its influence on SR yet. When d = 9 Å, the lower pressures of 50 and 100 MPa begin to increase the salt rejection. Lower
ð2Þ
Fig. 4. Influence of (a) pressure and (b) channel thickness on salt rejection of the GE nanochannel. Simulations in (b) is carried out under constant driving pressure DP = 300 MPa.
104
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
Fig. 5. Influence of electric field E (1–5 V/Å) on (a) velocity of water molecules inside the channel, (b) velocity of ions (average of Na+ and Cl) inside the channel, (c) number of ions inside the channel, as well as on salt rejection. Pressure difference is constant at DP = 300 MPa.
pressure under 50 MPa is beyond our scope in this study, but it’s worth investigating in future work based on the observation here. Different from the channel formed by GE multilayers, the salt rejection of porous GE single layer depends greatly on the pressure difference in a similar range from 100 to 225 MPa [18]. As for the influence of channel thickness, the most obvious observation from Fig. 4 is the sudden decrease from d = 8 Å to d = 9 Å, where water molecules inside the nanochannel change from one-layer to twolayer structure. However, when the channel thickness changes from 11 Å to 12 Å and the layer of water molecules transforms into three, there is only a small decrease. Therefore, d = 8 Å seems to be the critical d that can block small ions such as Na+ and Cl, but it should be noted that the water flux is also small at this d. In real application, there needs to be a balance of water flux and salt rejection for the design of GE derived membranes. For example, if the salt rejection at d = 9–13 Å still varies little under the applied pressure, it would be better to choose larger d for larger flux. At DP = 200 MPa, the critical diameter to block the salt ion is 5.5 Å for porous GE single layer (SR = 90% by Cohen-Tanugi and Grossman) [18], and 10.9 Å for CNT (SR = 30–60% by Thomas and Corry) [44].
exceeds 3 V/Å, while at d = 20 Å, SR goes beyond 90% only when E exceeds 6 V/Å. At d = 13 Å, there is a rapid increase of SR from 56% to 94% when E increases from 2 to 3 V/Å. Similar result is found at d = 10–12 Å, with SR increasing rapidly from 57 to 75% to over 90%. These observations could be explained by investigating the potential of mean force inside the channel and friction force between the ions and GE layers, which will be discussed in our future work. Electric field drives ions towards the GE channel walls, as shown in Fig. 6(a) and (c), and ions lose part of the hydration shell, same as those inside CNT within electric field [20]. We plot radius distribution function (RDF) profiles of chlorine-oxygen (Cl-O), sodium-oxygen (Na+-O) and oxygen-oxygen (O-O) pairs in Fig. 7.
3.3. Influence of external electric field When the electric field with intensity E ranging from 1 to 5 V/Å is added vertically to the GE nanochannel, the velocities of both water molecules and ions decrease drastically, as shown in Fig. 5. From 1 to 2 V/Å, the velocities decrease rapidly for both water ðIONÞ can nearly molecules and ions. It can also be observed that v reach zero at certain E depending on the channel thickness. Besides, more ions can be held inside the GE nanochannel with a larger electric field intensity, indicating that ions can be trapped inside the channel. Correspondingly, salt rejection can yield higher values, and SR can reach nearly 100% at high E. Also, the channel thickness plays an important role on the trend of SR change. For example, at d = 9 and 13 Å, SR reaches a value over 90% when E
Fig. 6. Snapshots of ions inside bilayer GE nanochannel when (a) channel thickness d = 13 Å, driving pressure DP = 300 MPa, and electric field E = 3 V/Å; (b) d = 13 Å, DP = 10 MPa, and E = 0.1 V/Å; (c) d = 20 Å, DP = 300 MPa, and E = 3 V/Å; (d) d = 20 Å, DP = 10 MPa, and E = 0.1 V/Å.
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
The first peak of RDF for O-O and Na+-O occurs at r = 2.82 and 2.46 Å respectively, consistent with previous research [25,45–47]. However, the first peak of Cl-O occurs at r = 3.42 Å, slightly larger than the reported range 3.13–3.35 Å [48–50]. When no electric field is imposed, the first RDF peak at d = 13 Å is always higher than that at d = 20 Å, indicating the GE nanochannel itself influences the distribution of water molecules as well as hydration structures around ions inside the channel. The RDF curves also vary greatly with electric field strength, especially for those of the ions. Taking pairs of Na+-O and Cl-O as examples, a second peak occurs apparently at d = 13 Å, when there is no electric field, and the peak decreases when the electric field is added. Therefore, the structure of hydration shell around ions is under the influence of both nanochannel and electric field. The distances of Na+ to the upper GE channel wall, Cl to the lower GE channel wall, as well as oxygen atom to the upper and lower GE channel walls are also examined. Fig. 8 highlights that the distance of ions to the GE layers decreases as the electric field is strengthened, but the distance stays the same for oxygen atoms. Therefore, the interaction between the ions and GE bilayers gets stronger, and the friction force to the flow increases. Because Cl has a larger radius than Na+, its distance to the GE layer is also larger. The results above also indicate that the applied electric field can separate ions from the feed flow and trap them inside the channel
105
to produce desalinated permeate. The model resembles the electrodes of downscaling capacitive deionization (CDI). Although we did not construct the model of porous GE layers, the results above can still show that layered GE electrodes have the potential for CDI applications. However, in practical applications use or laboratory experiments, the feed water velocity cannot be as large as 25 m/s which is the velocity of water molecules at d = 9 and 13 Å, DP = 300 MPa, and E = 3 V/Å, and the electric field cannot be this strong, as mentioned in Section 2. Therefore, cases at lower pressure (DP = 5 and 10 MPa) and lower electric field intensity (E = 0.01–0.1 and 0.05– 0.2 V/Å) at d = 13 and 20 Å are also tested. The results in Fig. 9 show different patterns compared with the results in Fig. 5. Velocities fall to several m/s and numbers of ions are also less than those at the high-pressure range with the numbers of water molecules stay the same. The most obvious difference is that the numbers of ions decrease as the electric field gets stronger, contrary with the results in Fig. 5. Meanwhile, the increase of salt rejection results from the decreased ion numbers rather than the velocity difference between ions and water molecules. This means that instead of trapping the ions inside the channel, the GE nanochannel with electric field repulses ions outside the channel instead. The ions inside the nanochannel still move forward along the channel wall as shown in Fig. 6(b) and (d), and ions stay farther away from the GE channel wall than those under a stronger electric field.
Fig. 7. RDF profiles of O-O, Na+-O and Cl-O inside the channel at channel thickness d = 13 and 20 Å, driving pressure DP = 300 MPa, and electric field intensity E = 0–5 V/Å, as well as bulk flow outside the channel.
106
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
Fig. 8. Average distance of Na+ to upper GE layer, Cl to lower GE layer, oxygen atom to upper GE layer as well as oxygen atom to lower GE layer inside the channel, when channel thickness d = 13 and 20 Å, pressure difference DP = 300 MPa, and electric field intensity E = 0–5 V/Å.
Fig. 9. (a–c) Influence of electric field E (0.05–0.2 V/Å) on number of ions (average of Na+ and Cl), velocity of water molecules and ions (average of Na+ and Cl) inside the channel, and salt rejection at constant pressure difference DP = 10 MPa. (d–f) Similar results as in a–c but for the electric field of 0.01–0.1 V/Å and the constant pressure difference DP = 5 MPa.
From the results in Fig. 9(c) and (f), it can be observed that after applying an electrical field of 0.05–0.2 V/Å and pressure of 10 MPa, the salt rejection rate for the GE bilayer channel with the widths of 20 Å and 13 Å increase by 40% and 100% respectively. The increase in the salt rejection rate can be explained by the results shown in Fig. 9(a) and (d). As the strength of the external electrical field increases, less ions enter into the channel, indicating that the physical filtration capability of the membrane is enhanced. In addition, the mobility of the ions inside the channel is also reduced, particularly for the large interlayer width of 20 Å. Therefore, the salt rejection rate increases correspondingly. Due to the lesser number of ions inside the channel, the potential of internal fouling of the membrane is also expected to be reduced. At the same time, it can be observed from Fig. 9(b) that the velocities of the water molecules inside the channel are also reduced for the large interlayer width of 20 Å at 10 MPa, indicating a reduced permeate flux. However, for the smaller interlayer width
of 13 Å or a weaker electrical field (Fig. 9b and e), the permeate flux is maintained. In summary, by applying the electric field across the GE-based membrane, better permeate quality and lesser ion channel fouling can be expected at the same driving pressure without sacrificing the permeate flux. Whether the energy consumption can be reduced together depends on how the electric field is imposed for application. For example, lower energy consumption can be expected if the electric field is introduced by functionalizing the GE membranes with the ionized oxygen functional groups [51]. 4. Conclusions In this study, molecular dynamics simulations are conducted to simulate the pressure-driven transport of salty water through the nanochannel formed by the GE bilayer. Simulations are carried out at the channel thickness ranging from 7 to 20 Å with and
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
without a vertical electric field. The number of water molecules and ions inside the nanochannel, as well as their velocities and distributions are investigated to understand their transport behaviors. The results show that the number of water molecules inside the nanochannel is independent of the driving pressure; while the number of ions inside the nanochannel is sensitive to the pressure change when the channel thickness d is small. Moreover, there is a more rapid increase of the number of ions than water molecules inside the nanochannel when d changes from 8 to 9 Å, where the structure of water molecules inside the channel changes from one-layer to two-layer. The velocities of ions are smaller than water molecules inside the channel, and such velocity difference reduces as d gets larger. Based on the different behavior of water molecules and ions inside the nanochannel, salt rejection is defined and calculated in terms of the velocity and number of water molecules and ions inside the channel. It is found that within the driving pressure range of 50–300 MPa, the pressure only influences the salt rejection at a small channel thickness of 9 Å when no electric field is imposed. Besides, at the high pressure of 300 MPa, the salt rejection decreases rapidly from 90% to 58% as d increases from 8 to 9 Å, and decreases slowly from 58% to 56% as d further increases from 9 to 13 Å. When an electric field is added vertically to the nanochannel with a large intensity ranging from 1 to 5 V/Å under the high driving pressure of 300 MPa, ions can be trapped inside the nanochannel when the electric field is sufficiently strong. Cases at lower pressures (5 and 10 MPa) and electric field intensities (0.01–0.1 and 0.05–0.25 V/Å) are also tested. It is found that the GE nanochannel repulses the ions outside the channel instead of trapping them inside. These results imply that improved permeate quality and lesser ion channel fouling can be expected without sacrificing the permeate flux when a vertical electric field is imposed. Overall, the simulation results in the present study can provide guidance for the design and application of GE-based RO membranes and CDI electrodes. References [1] C.Y. Chen, F. Yu, H.M. Zhou, J.H. Chen, J. Ma, Chem. J. Chin. Univ. – Chin. 36 (2015) 2516. [2] Y. Han, Z. Xu, C. Gao, Adv. Funct. Mater. 23 (2013) 3693. [3] H.B. Huang, Z.G. Song, N. Wei, L. Shi, Y.Y. Mao, Y.L. Ying, L.W. Sun, Z.P. Xu, X.S. Peng, Nat. Commun. 4 (2013) 1. [4] R.K. Joshi, P. Carbone, F.C. Wang, V.G. Kravets, Y. Su, I.V. Grigorieva, H.A. Wu, A. K. Geim, R.R. Nair, Science 343 (2014) 752. [5] R.R. Nair, H.A. Wu, P.N. Jayaram, I.V. Grigorieva, A.K. Geim, Science 335 (2012) 442. [6] S.C. O’Hern, M.S.H. Boutilier, J.C. Idrobo, Y. Song, J. Kong, T. Laoui, M. Atieh, R. Karnik, Nano Lett. 14 (2014) 1234. [7] S. Porada, L. Borchardt, M. Oschatz, M. Bryjak, J.S. Atchison, K.J. Keesman, S. Kaskel, P.M. Biesheuvel, V. Presser, Energy Environ. Sci. 6 (2013) 3700.
107
[8] X.T. Xu, Y. Liu, M. Wang, X.X. Yang, C. Zhu, T. Lu, R. Zhao, L.K. Pan, Electrochim. Acta 188 (2016) 406. [9] W. Choi, J. Choi, J. Bang, J.H. Lee, ACS Appl. Mater Interf. 5 (2013) 12510. [10] H.W. Kim, H.W. Yoon, S.M. Yoon, B.M. Yoo, B.K. Ahn, Y.H. Cho, H.J. Shin, H. Yang, U. Paik, S. Kwon, J.Y. Choi, H.B. Park, Science 342 (2013) 9. [11] N. Wei, X.S. Peng, Z.P. Xu, ACS Appl. Mater. Interf. 6 (2014) 5877. [12] B.X. Mi, Science 343 (2014) 740. [13] A. Aghigh, V. Alizadeh, H.Y. Wong, M.S. Islam, N. Amin, M. Zaman, Desalination 365 (2015) 389. [14] L. Yu, Y. Zhang, B. Zhang, J. Liu, H. Zhang, C. Song, J. Membr. Sci. 447 (2013) 452. [15] K. Celebi, J. Buchheim, R.M. Wyss, A. Droudian, P. Gasser, I. Shorubalko, J.I. Kye, C. Lee, H.G. Park, Science 344 (2014) 289. [16] A. Akbari, P. Sheath, S.T. Martin, D.B. Shinde, M. Shaibani, P.C. Banerjee, R. Tkacz, D. Bhattacharyya, M. Majumder, Nat. Commun. 7 (2016) 10891. [17] D. Cohen-Tanugi, L.C. Lin, J.C. Grossman, Nano Lett. 16 (2016) 1027. [18] D. Cohen-Tanugi, J.C. Grossman, Nano Lett. 12 (2012) 3602. [19] J.G. Gai, X.L. Gong, W.W. Wang, X. Zhang, W.L. Kang, J. Mater. Chem. A 2 (2014) 4023. [20] Y. Kang, Z.S. Zhang, H. Shi, J.Q. Zhang, L.J. Liang, Q. Wang, H. Agren, Y.Q. Tu, Nanoscale 6 (2014) 10666. [21] S. Rikhtehgaran, A. Lohrasebi, Desalination 365 (2015) 176. [22] K. Sint, B. Wang, P. Kral, J. Am. Chem. Soc. 130 (2008) 16448. [23] B. Liu, R.B. Wu, J.A. Baimova, H. Wu, A.W.K. Law, S.V. Dmitriev, K. Zhou, Phys. Chem. Chem. Phys. 18 (2016) 1886. [24] Q. Tu, Q. Yang, H. Wang, S. Li, Sci. Rep. 6 (2016) 26183. [25] A.T. Nasrabadi, M. Foroutan, Desalination 277 (2011) 236. [26] D. Zhang, T. Yan, L. Shi, Z. Peng, X. Wen, J. Zhang, J. Mater. Chem. 22 (2012) 14696. [27] M. Noked, A. Soffer, D. Aurbach, J. Solid State Electrochem. 15 (2011) 1563. [28] M.E. Suss, S. Porada, X. Sun, P.M. Biesheuvel, J. Yoon, V. Presser, Energy Environ Sci. 8 (2015) 2296. [29] S. Banerjee, S. Murad, I.K. Puri, Chem. Phys. Lett. 434 (2007) 292. [30] S. Plimpton, J. Comput. Phys. 117 (1995) 1. [31] D.J. Price, C.L. Brooks, J. Chem. Phys. 121 (2004) 10096. [32] A.D. MacKerell, D. Bashford, M. Bellott, R.L. Dunbrack, J.D. Evanseck, M.J. Field, S. Fischer, J. Gao, H. Guo, S. Ha, D. Joseph-McCarthy, L. Kuchnir, K. Kuczera, F.T. K. Lau, C. Mattos, S. Michnick, T. Ngo, D.T. Nguyen, B. Prodhom, W.E. Reiher, B. Roux, M. Schlenkrich, J.C. Smith, R. Stote, J. Straub, M. Watanabe, J. Wiorkiewicz-Kuczera, D. Yin, M. Karplus, J. Chem. Phys. B 102 (1998) 3586. [33] I.S. Joung, T.E. Cheatham III, J. Chem. Phys. B 112 (2008) 9020. [34] E.M. Stuve, Chem. Phys. Lett. 519 (2012) 1. [35] G. Cicero, J.C. Grossman, E. Schwegler, F. Gygi, G. Galli, J. Am. Chem. Soc. 130 (2008) 1871. [36] J.A. Thomas, A.J.H. McGaughey, Phys. Rev. Lett. 102 (2009) 4. [37] A.A. Chialvo, P.T. Cummings, J. Phys. Chem. A 115 (2011) 5918. [38] D. Konatham, J. Yu, T.A. Ho, A. Striolo, Langmuir 29 (2013) 11884. [39] D.G. Bounds, Mol. Phys. 54 (1985) 1335. [40] H. Ohtaki, T. Radnai, Chem. Rev. 93 (1993) 1157. [41] Y. Marcus, J. Chem. Soc. Faraday Trans. 87 (1991) 2995. [42] T.A. Hilder, D. Gordon, S.H. Chung, Small 5 (2009) 2183. [43] K. Leung, S.B. Rempe, J. Comput. Theor. Nanos. 6 (2009) 1948. [44] M. Thomas, B. Corry, Philos. Trans. R. Soc. A – Math. Phys. Eng. Sci. 374 (2016) 20. [45] M. Cavallari, C. Cavazzoni, M. Ferrario, Mol. Phys. 102 (2004) 959. [46] P. Gallo, D. Corradini, M. Rovere, Phys. Chem. Chem. Phys. 13 (2011) 19814. [47] Q. Shao, J. Zhou, L.H. Lu, X.H. Lu, Y.D. Zhu, S.Y. Jiang, Nano Lett. 9 (2009) 989. [48] F. Bruni, S. Imberti, R. Mancinelli, M.A. Ricci, J. Phys. Chem. 136 (2012) 1. [49] R. Caminiti, G. Licheri, G. Piccaluga, G. Pinna, Chem. Phys. Lett. 47 (1977) 275. [50] Z.H. Duan, Z.G. Zhang, Mol. Phys. 101 (2003) 1501. [51] P. Sun, R. Ma, W. Ma, J. Wu, K. Wang, T. Sasaki, H. Zhu, NPG Asia Mater. 8 (2016) 259.
Contents lists available at ScienceDirect
Computational Materials Science journal homepage: www.elsevier.com/locate/commatsci
Transport of salty water through graphene bilayer in an electric field: A molecular dynamics study Hui Zhang a, Bo Liu a, Mao-See Wu b, Kun Zhou a,b,⇑, Adrian Wing-Keung Law a,c,⇑ a Environmental Process Modelling Centre, Nanyang Environment & Water Research Institute (NEWRI), Nanyang Technological University, 1 Cleantech Loop, CleanTech One, 637141, Singapore b School of Mechanical and Aerospace Engineering, Nanyang Technological University, 50 Nanyang Avenue, Singapore 639798, Singapore c School of Civil and Environmental Engineering, Nanyang Technological University, 50 Nanyang Avenue, Singapore 639798, Singapore
a r t i c l e
i n f o
Article history: Received 1 November 2016 Received in revised form 24 January 2017 Accepted 28 January 2017
Keywords: Molecular dynamics Graphene bilayer Electric field Salt rejection
a b s t r a c t Molecular dynamics simulations are carried out to study the pressure-driven transport of salty (NaCl) water through the nanochannels formed by a graphene (GE) bilayer with and without a vertical electric field. The simulation results show that the channel thickness influences not only the velocity and numbers of water molecules and ions inside the channel, but also their responses to varied pressure and electric field intensity. The salt rejection rate is characterized by considering both the numbers of ions and water molecules inside the channel and their velocities. When the electric field is not imposed, the salt rejection rate varies slightly in the high-pressure range (50–300 MPa), and is mainly affected by the channel thickness. When an electric field with high intensity (1–5 V/Å) is imposed vertically under high driving pressure (300 MPa), ions can be trapped inside the nanochannel. Simulations under lower electric field (0.01–0.1 V/Å and 0.05–0.2 V/Å) and pressure (5 MPa and 10 MPa) are also conducted. Under these circumstances, GE nanochannels can reject ions from entering the channel, indicating improved salt rejection rate and reduced ion channel accumulation. Ó 2017 Elsevier B.V. All rights reserved.
1. Introduction Graphene (GE) has recently attracted substantial attention in water desalination techniques [1–8]. Reverse osmosis (RO) membranes derived from GE have been reported to show advantages over conventional RO membranes because of ultimate thinness, high mechanical strength, and more importantly, high water permeability [2–5,9–11]. GE derived membranes can be fabricated by stacking functionalized GE nanosheets, where nanochannel and nanopores are formed. Interlayer distance of graphene oxide (GO) membrane can be less than 0.7 nm if made by vacuum filtration in the dry state, enabling water molecules to permeate through while blocking hydrated ions and other target molecules [5,3,12]. Confined in the nanostructures, water still shows high permeability through GO membranes. This is attributed to the capillary network formed by the pristine GE domains, where water can flow almost frictionlessly [5,11].
⇑ Corresponding authors at: Environmental Process Modelling Centre, Nanyang Environment & Water Research Institute (NEWRI), Nanyang Technological University, 1 Cleantech Loop, CleanTech One, 637141, Singapore. E-mail addresses: [email protected] (K. Zhou), [email protected] (A.W.-K. Law). http://dx.doi.org/10.1016/j.commatsci.2017.01.039 0927-0256/Ó 2017 Elsevier B.V. All rights reserved.
GE-based membranes also have superior energy efficiency compared with other RO technology [13]. Meanwhile, functionalized GO can enhance the hydrophilicity of the membrane and provide negative charges, and thus possessing anti-fouling properties [14]. However, the scale-up ability remains a critical issue for the mass production of GE-based membranes. Many studies have reported on various possible fabrication approaches, such as chemical vapour deposition and subsequent transfer to substrates [15], restacking GO flakes by filtration on a back support [2], etc. Most of these methods suffer from scalability issues. Recently, Akbari et al. developed a method to produce large-area, highly ordered, and multi-layered GE-based films by making use of the shear alignment of discotic nematic GO fluid [16]. This scalable method brings hope to the mass production of GE-based membranes for industrial applications. Molecular dynamics (MD) simulation has been used to investigate the transport of salty water through GE nanopores and nanochannels [4,11,17–24]. GE single layer with nanopores bonded with chemical functional groups at the pore ring can limit the salt passage depending on the pore size, functional groups and applied pressure or electric field [18–22]. Cohen-Tanugi et al. also carried out MD simulations of salty water passing through multilayer nanoporous GE with different interlayer width as well
101
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
as pore offset, and demonstrated that salt separation performance can be manipulated by different multilayer design [17]. Simulation results also showed that the water molecules display layered structure inside GE nanochannels which affects the water flow rate significantly [4,8]. GE nanochannels with the thickness of allowing one monolayer of water molecules passing through do not allow ions inside [4]. However, a systematic study of ion transport through the nanochannel formed by GE bilayers, especially the interaction of ions with water molecules and GE, is still needed to provide further guidance on the design and application of GEbased RO membranes. High surface area and excellent conductivity also make GE feasible for the fabrication of electrodes in capacitive deionization (CDI) [1,8,26,27]. The experimental results from these studies showed that the electrosorption capacity of GE-derived electrodes can be enhanced by mixing GE with other chemicals or creating specific structures. In CDI, salty water flows between the electrodes [26,27]. MD simulations can provide insights into the effect of nanostructure on CDI performance which is still not fully determined or optimized [28]. The present work focuses on the behavior of ions and water molecules in flowing solutions rather than static solutions, while studies so far often use the latter one for MD simulations due to their relative ease to be set up [25,29]. Nasrabadi and Foroutan carried out MD simulations to demonstrate that carbon nanotube (CNT) electrodes are suitable for ion separation by investigating the atomic interaction among ions, water molecules and CNT [25]. However, the feasibility study of GE derived electrodes is mostly based on material characterization and experiments so far. The present study is targeted at a systematic investigation of the transport of salty water through the nanochannel formed by a GE bilayer with different channel thickness, applied pressure as well as external electric field using MD simulations. Interactions among the ions, water molecules and GE are also investigated. The simulation results provide useful guidelines for the design and application of GE-derived material for RO or CDI desalination.
2. Simulation details The studied system consists of a nanochannel formed by the bilayer GE and two water reservoirs at each end of the channel constrained by GE walls, as shown in Fig. 1. The left-hand feeding
reservoir is initially filled with salty water containing 16 Na+ and 16 Cl dissolved with 922 water molecules, and the right-hand receiving reservoir with pure water containing 954 water molecules. The salinity in the feeding reservoir is 4.86 g/L, slightly higher than sea water salinity (3.5 g/L). Both reservoirs have the same dimension of 40 20 40 Å3. The GE bilayer is parallel to the xy plane with a constant length L 32 Å, while its thickness d varies from 7 to 20 Å. The GE constrain walls are set normal to the bilayer plane to restrict water and ion movement. Simulations are carried out using the large-scale atomic/molecular massively parallel simulator (LAMMPS) package [30]. Boundary conditions are set as periodic in all three dimensions to maintain the continuous feeding. The Transferable Intermolecular Potential 3P (TIT3P) [31] is applied to model the water molecules, and partial charges of O and H are 0.834e and 0.417e. Other atomic interactions are modelled using CHARMM force field [32], including Coulombics and van der Waals (vdW) interactions. The long-range Coulombics are computed with the particle-particleparticle-mesh (PPPM) method with a cutoff of 12 Å, and vdW interaction is modelled using Lennard-Jones (LJ) potential with a cutoff at 10 Å. LJ parameters and charges of each atom species are summarized in Table 1 according to previous research studies [23,32,33]. The LJ pair parameters between species i and j are calculated using the Lorentz-Berthelot combining rules except for C-O and C-H interactions. The energy well e and zero-across distance r between C and O are set as 0.0937 kcal/mol and 3.19 Å, while e and r between C and H are all set as 0 [23]. The NVT ensemble with a Nosé-hoover thermostat is used to keep the temperature at 300 K. The GE nanochannel and constrain walls are treated as rigid body with carbon atoms fixed. The SHAKE algorithm is used to constrain the bond and angle of water molecules to maintain the rigid molecule structure. After all the water molecules are relaxed for 50 ps with a timestep of 0.5 fs to reach equilibrium state, a pressure difference DP is applied by adding a force f on all the O atoms in the left feeding reservoir along the positive X-direction. The force f is calculated as [23]
f ¼ DPLy Lz =n
ð1Þ
where n is the number of water molecules in the left water reservoir, and LyLz is the cross-section area shown in Fig. 1. Our simulations focus on the pressure difference DP ranging from 50 to 300 MPa to obtain data for statistical calculations. The external electric field is added in the channel region and vertical to the xy plane. Because ions and water molecules move very fast due to the high pressure, the applied electric field E also has a large intensity (1–5 V/Å) to obtain distinguishable results. The intensity range here is with the same magnitude of electric field that can break bonds and induce ionization of water molecules [34], but the rigid structure of water molecules is kept using the SHAKE algorithm. Cases with lower pressures (5 MPa and 10 MPa) and lower electric field intensities (0.01–0.1 V/Å and 0.05–0.20 V/Å) are also tested, which are discussed in Section 3.3. The simulated system can reach the steady state within around 100 ps with the imposed force f and electric field E. The simulation is then carried on for another 1.5 ns (2.5 ns when DP = 5 and 10 MPa) to collect data. Table 1 LJ parameters and charges used in the simulation.
Fig. 1. Front and top view of the simulated system consisting of two reservoirs and a nanochannel formed by the GE bilayer. At the initial condition, the left-hand reservoir is filled with salt (Na+ and Cl) water, and the right-hand reservoir is filled with pure water. The two reservoirs are separated by GE nanochannel and constrain wall.
Species
e (kcal/mol)
r (Å)
q (e)
H O C Na+ Cl
0.0461 0.1521 0.0700 0.0874 0.0356
0.40 3.15 3.98 2.74 5.03
0.417 0.834 0 1 1
102
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
3. Results and discussion 3.1. Influence of pressure and channel thickness of In each simulation, the number N and average velocity v water molecules and ions are obtained after the system reaches remain constant with slight fluctuthe steady state where N and v of water molecules and ions are compared ations. The N and v under the applied pressure difference DP varying from 50 to 300 MPa at the channel thickness d of 9 Å, 13 Å and 20 Å, as shown in Fig. 2. The three thicknesses chosen here correspond with twolayer, three-layer and bulk water distribution structures inside the nanochannel, reported by previous studies [4,23,35]. The one-layer water flow at smaller thickness is not chosen because ions can hardly enter such narrow channel [4,17]. An almost linear increase of velocity in response to pressure DP can be observed for both water molecules and ions in Fig. 2 when the electric field is not imposed. At the same channel thickness, the velocities of ions have similar magnitude to but are slightly smaller than the velocities of water molecules. The slight difference of velocities between the water molecules and ions becomes more apparent when the channel thickness is smaller. At a larger channel thickness, the velocities of both water molecules and ions increase more rapidly with the increasing pressure, suggesting that a weaker constraint by the nanochannel enables more significant influence of the driving pressure. The velocity of water molecules inside the CNTs with diameters of 11–16.6 Å at DP = 75 MPa varies from 2 to 6.5 m/s [36]. In the present simulations, the velocities of water inside the GE nanochannel of d = 13 Å are 10.6 m/s and 16.6 m/s at DP = 50 MPa and DP = 100 MPa respectively, higher than that in the CNTs. It can be observed from Fig. 2 that the number of water molecules inside the channel varies little in the pressure range of 50– 300 MPa. As for the number of ions, at d = 20 Å, no obvious change
occurs. At d = 13 Å, the lower pressure of 50 MPa starts to show its effect on reducing N (ION); at d = 9 Å, a more obvious drop of N (ION) can be observed. This suggests that the number of ions inside the nanochannel is more sensitive to pressure changes than water molecules, especially at smaller channel thickness. Dissolved in water, ions are surrounded with hydration shell and they tend to lose part of their hydration shell when confined in other nanostructures [37,38]. In an aqueous NaCl solution of similar salinity with that in this study, the hydration radii of Cl and Na+ are 3.23 Å and 2.35 Å respectively [39,40]. The equilibrium vdW distance between water molecules on the hydration shell and the GE layers is 3.19 Å (equal to rC-O), as illustrated in Section 2. Therefore, in order for ions to get into the narrow nanochannel at d = 9 Å, their hydration shell has to be partly peeled off, which demands more energy. The dependence of transport properties of water and ions inside the channel on its thickness d is described in Fig. 3 with a constant pressure difference and without electric field. Our previous work has reported the dependence of N (H2O) on the channel thickness for d = 7–13 Å in detail, including the sharp increase between d = 8 Å and 9 Å as well as the less sharp increase between d = 11 Å and 12 Å [23]. In this study, d = 20 Å is also investigated to indicate the dependence of d over a wider range. For d = 7 Å, ions cannot pass through. This can be explained by the radii of Na+ and Cl (1.02 Å and 1.81 Å) [41] and their equilibrium vdW distance with the GE layer (rNa-C = 3.36 Å and rCl-C = 4.50 Å). The increase in gradient of N (ION) is similar to that of N (H2O), except there is a sharper increase from d = 8 Å to d = 9 Å for N (ION), where water molecules change from two-layer to threelayer. From the velocity profile, it can be observed that although v (H2O) and v (ION) increase with the similar trend as d gets larger, v (H2O) is always larger than v (ION). The velocity difference indicates that the constraining effect of the nanochannel on ions is stronger than that on water molecules. The difference between
Fig. 2. Influence of pressure on (a) velocity of water molecules, (b) ratio of velocity of ions (average of Na+ and Cl) vs velocity of water molecules, (c) number of water molecules, and (d) number of ions (sum of Na+ and Cl) inside the GE nanochannel of thickness d = 9, 13 and 20 Å when electric field is not applied.
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
103
Fig. 3. Influence of channel thickness on (a) number of water molecules inside the channel, (b) number of ions inside the channel, and (c) velocity of water molecules and ions (average of Na+ and Cl) inside the channel at constant pressure difference DP = 300 MPa.
v (H2O) and v (ION) reduces as d gets larger, because a nanochannel with larger d has a weaker constraining effect.
channel. Therefore, in the present study, the salt rejection SR is defined as below,
3.2. Salt rejection
SR ¼ ð1 KC=C 0 Þ 100%
Salt rejection rate (SR) is an important parameter to assess the performance of RO membrane. In engineering, SR is defined as one minus the ratio of permeate salinity versus feed salinity. In MD simulations of either GE or CNT, the definition of SR generally follows the engineering approach. Previous research calculates SR by comparing the number of permeated water molecules and ions to the initial values in the feed reservoir, or only by counting the number of permeated ions, regardless of whether that boundary condition is periodic or not [17,42–44]. Periodic boundary condition enables the water molecules and ions to move from one reservoir to another without passing through the membrane. Therefore, it brings limitations to the general definition of SR based on counting numbers, which is also highlighted by Thomas and Corry [44]. However, under non-periodic boundary condition, another problem arises regarding when to count the numbers. This is because ion concentrations under non-periodic boundary conditions change from time to time rather than reaching a steady state. Thus, some researches use the time point instead when certain amount of water molecules or ions enter into the permeate reservoir [17,43]. As discussed in Section 3.1, the GE nanochannel can influence the ions passing through by reducing their velocity inside the channel as well as the opportunity for them to enter the channel. These two aspects can be expressed by the comparison of velocity and number between the ions and the water molecules inside the
where K ¼ v ðIONÞ=v ðH2 OÞ is the ratio of average velocities of ions and water molecules inside the channel; C = N (ION)/N (H2O), the ion concentration inside the channel, expressed as the ratio of the average numbers of ions and water molecules inside the channel; and C 0 ¼ N0 ðIONÞ=N 0 ðH2OÞ, the initial ion concentration expressed as the ratio of the initial numbers of ions and water molecules in the left feeding reservoir. The definition of SR here can overcome the difficulties of counting numbers under the periodic boundary condition and help decide on the time point under non-periodic boundary condition. Each simulation first runs for a rather long time so that the system reaches the steady state, simulation data are then collected to calculate the average values of each parameter in Eq. (2) to make SR more statistically reasonable. In fact, given any time period, Eq. (2) expresses the same meaning as that of the general definition by integration over time. In this manner, we can ensure that Eq. (2) calculates a steady-state SR. Indeed, we are not able to obtain large salinity difference at each side of the channel, but under high pressure range, its influence on diffusion could be neglected. Fig. 4 shows the salt rejection under the influence of pressure difference DP and channel thickness d. Within the investigated range of DP, SR at d = 13 and 20 Å under different pressure does not vary significantly, implying that pressure has not started to show its influence on SR yet. When d = 9 Å, the lower pressures of 50 and 100 MPa begin to increase the salt rejection. Lower
ð2Þ
Fig. 4. Influence of (a) pressure and (b) channel thickness on salt rejection of the GE nanochannel. Simulations in (b) is carried out under constant driving pressure DP = 300 MPa.
104
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
Fig. 5. Influence of electric field E (1–5 V/Å) on (a) velocity of water molecules inside the channel, (b) velocity of ions (average of Na+ and Cl) inside the channel, (c) number of ions inside the channel, as well as on salt rejection. Pressure difference is constant at DP = 300 MPa.
pressure under 50 MPa is beyond our scope in this study, but it’s worth investigating in future work based on the observation here. Different from the channel formed by GE multilayers, the salt rejection of porous GE single layer depends greatly on the pressure difference in a similar range from 100 to 225 MPa [18]. As for the influence of channel thickness, the most obvious observation from Fig. 4 is the sudden decrease from d = 8 Å to d = 9 Å, where water molecules inside the nanochannel change from one-layer to twolayer structure. However, when the channel thickness changes from 11 Å to 12 Å and the layer of water molecules transforms into three, there is only a small decrease. Therefore, d = 8 Å seems to be the critical d that can block small ions such as Na+ and Cl, but it should be noted that the water flux is also small at this d. In real application, there needs to be a balance of water flux and salt rejection for the design of GE derived membranes. For example, if the salt rejection at d = 9–13 Å still varies little under the applied pressure, it would be better to choose larger d for larger flux. At DP = 200 MPa, the critical diameter to block the salt ion is 5.5 Å for porous GE single layer (SR = 90% by Cohen-Tanugi and Grossman) [18], and 10.9 Å for CNT (SR = 30–60% by Thomas and Corry) [44].
exceeds 3 V/Å, while at d = 20 Å, SR goes beyond 90% only when E exceeds 6 V/Å. At d = 13 Å, there is a rapid increase of SR from 56% to 94% when E increases from 2 to 3 V/Å. Similar result is found at d = 10–12 Å, with SR increasing rapidly from 57 to 75% to over 90%. These observations could be explained by investigating the potential of mean force inside the channel and friction force between the ions and GE layers, which will be discussed in our future work. Electric field drives ions towards the GE channel walls, as shown in Fig. 6(a) and (c), and ions lose part of the hydration shell, same as those inside CNT within electric field [20]. We plot radius distribution function (RDF) profiles of chlorine-oxygen (Cl-O), sodium-oxygen (Na+-O) and oxygen-oxygen (O-O) pairs in Fig. 7.
3.3. Influence of external electric field When the electric field with intensity E ranging from 1 to 5 V/Å is added vertically to the GE nanochannel, the velocities of both water molecules and ions decrease drastically, as shown in Fig. 5. From 1 to 2 V/Å, the velocities decrease rapidly for both water ðIONÞ can nearly molecules and ions. It can also be observed that v reach zero at certain E depending on the channel thickness. Besides, more ions can be held inside the GE nanochannel with a larger electric field intensity, indicating that ions can be trapped inside the channel. Correspondingly, salt rejection can yield higher values, and SR can reach nearly 100% at high E. Also, the channel thickness plays an important role on the trend of SR change. For example, at d = 9 and 13 Å, SR reaches a value over 90% when E
Fig. 6. Snapshots of ions inside bilayer GE nanochannel when (a) channel thickness d = 13 Å, driving pressure DP = 300 MPa, and electric field E = 3 V/Å; (b) d = 13 Å, DP = 10 MPa, and E = 0.1 V/Å; (c) d = 20 Å, DP = 300 MPa, and E = 3 V/Å; (d) d = 20 Å, DP = 10 MPa, and E = 0.1 V/Å.
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
The first peak of RDF for O-O and Na+-O occurs at r = 2.82 and 2.46 Å respectively, consistent with previous research [25,45–47]. However, the first peak of Cl-O occurs at r = 3.42 Å, slightly larger than the reported range 3.13–3.35 Å [48–50]. When no electric field is imposed, the first RDF peak at d = 13 Å is always higher than that at d = 20 Å, indicating the GE nanochannel itself influences the distribution of water molecules as well as hydration structures around ions inside the channel. The RDF curves also vary greatly with electric field strength, especially for those of the ions. Taking pairs of Na+-O and Cl-O as examples, a second peak occurs apparently at d = 13 Å, when there is no electric field, and the peak decreases when the electric field is added. Therefore, the structure of hydration shell around ions is under the influence of both nanochannel and electric field. The distances of Na+ to the upper GE channel wall, Cl to the lower GE channel wall, as well as oxygen atom to the upper and lower GE channel walls are also examined. Fig. 8 highlights that the distance of ions to the GE layers decreases as the electric field is strengthened, but the distance stays the same for oxygen atoms. Therefore, the interaction between the ions and GE bilayers gets stronger, and the friction force to the flow increases. Because Cl has a larger radius than Na+, its distance to the GE layer is also larger. The results above also indicate that the applied electric field can separate ions from the feed flow and trap them inside the channel
105
to produce desalinated permeate. The model resembles the electrodes of downscaling capacitive deionization (CDI). Although we did not construct the model of porous GE layers, the results above can still show that layered GE electrodes have the potential for CDI applications. However, in practical applications use or laboratory experiments, the feed water velocity cannot be as large as 25 m/s which is the velocity of water molecules at d = 9 and 13 Å, DP = 300 MPa, and E = 3 V/Å, and the electric field cannot be this strong, as mentioned in Section 2. Therefore, cases at lower pressure (DP = 5 and 10 MPa) and lower electric field intensity (E = 0.01–0.1 and 0.05– 0.2 V/Å) at d = 13 and 20 Å are also tested. The results in Fig. 9 show different patterns compared with the results in Fig. 5. Velocities fall to several m/s and numbers of ions are also less than those at the high-pressure range with the numbers of water molecules stay the same. The most obvious difference is that the numbers of ions decrease as the electric field gets stronger, contrary with the results in Fig. 5. Meanwhile, the increase of salt rejection results from the decreased ion numbers rather than the velocity difference between ions and water molecules. This means that instead of trapping the ions inside the channel, the GE nanochannel with electric field repulses ions outside the channel instead. The ions inside the nanochannel still move forward along the channel wall as shown in Fig. 6(b) and (d), and ions stay farther away from the GE channel wall than those under a stronger electric field.
Fig. 7. RDF profiles of O-O, Na+-O and Cl-O inside the channel at channel thickness d = 13 and 20 Å, driving pressure DP = 300 MPa, and electric field intensity E = 0–5 V/Å, as well as bulk flow outside the channel.
106
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
Fig. 8. Average distance of Na+ to upper GE layer, Cl to lower GE layer, oxygen atom to upper GE layer as well as oxygen atom to lower GE layer inside the channel, when channel thickness d = 13 and 20 Å, pressure difference DP = 300 MPa, and electric field intensity E = 0–5 V/Å.
Fig. 9. (a–c) Influence of electric field E (0.05–0.2 V/Å) on number of ions (average of Na+ and Cl), velocity of water molecules and ions (average of Na+ and Cl) inside the channel, and salt rejection at constant pressure difference DP = 10 MPa. (d–f) Similar results as in a–c but for the electric field of 0.01–0.1 V/Å and the constant pressure difference DP = 5 MPa.
From the results in Fig. 9(c) and (f), it can be observed that after applying an electrical field of 0.05–0.2 V/Å and pressure of 10 MPa, the salt rejection rate for the GE bilayer channel with the widths of 20 Å and 13 Å increase by 40% and 100% respectively. The increase in the salt rejection rate can be explained by the results shown in Fig. 9(a) and (d). As the strength of the external electrical field increases, less ions enter into the channel, indicating that the physical filtration capability of the membrane is enhanced. In addition, the mobility of the ions inside the channel is also reduced, particularly for the large interlayer width of 20 Å. Therefore, the salt rejection rate increases correspondingly. Due to the lesser number of ions inside the channel, the potential of internal fouling of the membrane is also expected to be reduced. At the same time, it can be observed from Fig. 9(b) that the velocities of the water molecules inside the channel are also reduced for the large interlayer width of 20 Å at 10 MPa, indicating a reduced permeate flux. However, for the smaller interlayer width
of 13 Å or a weaker electrical field (Fig. 9b and e), the permeate flux is maintained. In summary, by applying the electric field across the GE-based membrane, better permeate quality and lesser ion channel fouling can be expected at the same driving pressure without sacrificing the permeate flux. Whether the energy consumption can be reduced together depends on how the electric field is imposed for application. For example, lower energy consumption can be expected if the electric field is introduced by functionalizing the GE membranes with the ionized oxygen functional groups [51]. 4. Conclusions In this study, molecular dynamics simulations are conducted to simulate the pressure-driven transport of salty water through the nanochannel formed by the GE bilayer. Simulations are carried out at the channel thickness ranging from 7 to 20 Å with and
H. Zhang et al. / Computational Materials Science 131 (2017) 100–107
without a vertical electric field. The number of water molecules and ions inside the nanochannel, as well as their velocities and distributions are investigated to understand their transport behaviors. The results show that the number of water molecules inside the nanochannel is independent of the driving pressure; while the number of ions inside the nanochannel is sensitive to the pressure change when the channel thickness d is small. Moreover, there is a more rapid increase of the number of ions than water molecules inside the nanochannel when d changes from 8 to 9 Å, where the structure of water molecules inside the channel changes from one-layer to two-layer. The velocities of ions are smaller than water molecules inside the channel, and such velocity difference reduces as d gets larger. Based on the different behavior of water molecules and ions inside the nanochannel, salt rejection is defined and calculated in terms of the velocity and number of water molecules and ions inside the channel. It is found that within the driving pressure range of 50–300 MPa, the pressure only influences the salt rejection at a small channel thickness of 9 Å when no electric field is imposed. Besides, at the high pressure of 300 MPa, the salt rejection decreases rapidly from 90% to 58% as d increases from 8 to 9 Å, and decreases slowly from 58% to 56% as d further increases from 9 to 13 Å. When an electric field is added vertically to the nanochannel with a large intensity ranging from 1 to 5 V/Å under the high driving pressure of 300 MPa, ions can be trapped inside the nanochannel when the electric field is sufficiently strong. Cases at lower pressures (5 and 10 MPa) and electric field intensities (0.01–0.1 and 0.05–0.25 V/Å) are also tested. It is found that the GE nanochannel repulses the ions outside the channel instead of trapping them inside. These results imply that improved permeate quality and lesser ion channel fouling can be expected without sacrificing the permeate flux when a vertical electric field is imposed. Overall, the simulation results in the present study can provide guidance for the design and application of GE-based RO membranes and CDI electrodes. References [1] C.Y. Chen, F. Yu, H.M. Zhou, J.H. Chen, J. Ma, Chem. J. Chin. Univ. – Chin. 36 (2015) 2516. [2] Y. Han, Z. Xu, C. Gao, Adv. Funct. Mater. 23 (2013) 3693. [3] H.B. Huang, Z.G. Song, N. Wei, L. Shi, Y.Y. Mao, Y.L. Ying, L.W. Sun, Z.P. Xu, X.S. Peng, Nat. Commun. 4 (2013) 1. [4] R.K. Joshi, P. Carbone, F.C. Wang, V.G. Kravets, Y. Su, I.V. Grigorieva, H.A. Wu, A. K. Geim, R.R. Nair, Science 343 (2014) 752. [5] R.R. Nair, H.A. Wu, P.N. Jayaram, I.V. Grigorieva, A.K. Geim, Science 335 (2012) 442. [6] S.C. O’Hern, M.S.H. Boutilier, J.C. Idrobo, Y. Song, J. Kong, T. Laoui, M. Atieh, R. Karnik, Nano Lett. 14 (2014) 1234. [7] S. Porada, L. Borchardt, M. Oschatz, M. Bryjak, J.S. Atchison, K.J. Keesman, S. Kaskel, P.M. Biesheuvel, V. Presser, Energy Environ. Sci. 6 (2013) 3700.
107
[8] X.T. Xu, Y. Liu, M. Wang, X.X. Yang, C. Zhu, T. Lu, R. Zhao, L.K. Pan, Electrochim. Acta 188 (2016) 406. [9] W. Choi, J. Choi, J. Bang, J.H. Lee, ACS Appl. Mater Interf. 5 (2013) 12510. [10] H.W. Kim, H.W. Yoon, S.M. Yoon, B.M. Yoo, B.K. Ahn, Y.H. Cho, H.J. Shin, H. Yang, U. Paik, S. Kwon, J.Y. Choi, H.B. Park, Science 342 (2013) 9. [11] N. Wei, X.S. Peng, Z.P. Xu, ACS Appl. Mater. Interf. 6 (2014) 5877. [12] B.X. Mi, Science 343 (2014) 740. [13] A. Aghigh, V. Alizadeh, H.Y. Wong, M.S. Islam, N. Amin, M. Zaman, Desalination 365 (2015) 389. [14] L. Yu, Y. Zhang, B. Zhang, J. Liu, H. Zhang, C. Song, J. Membr. Sci. 447 (2013) 452. [15] K. Celebi, J. Buchheim, R.M. Wyss, A. Droudian, P. Gasser, I. Shorubalko, J.I. Kye, C. Lee, H.G. Park, Science 344 (2014) 289. [16] A. Akbari, P. Sheath, S.T. Martin, D.B. Shinde, M. Shaibani, P.C. Banerjee, R. Tkacz, D. Bhattacharyya, M. Majumder, Nat. Commun. 7 (2016) 10891. [17] D. Cohen-Tanugi, L.C. Lin, J.C. Grossman, Nano Lett. 16 (2016) 1027. [18] D. Cohen-Tanugi, J.C. Grossman, Nano Lett. 12 (2012) 3602. [19] J.G. Gai, X.L. Gong, W.W. Wang, X. Zhang, W.L. Kang, J. Mater. Chem. A 2 (2014) 4023. [20] Y. Kang, Z.S. Zhang, H. Shi, J.Q. Zhang, L.J. Liang, Q. Wang, H. Agren, Y.Q. Tu, Nanoscale 6 (2014) 10666. [21] S. Rikhtehgaran, A. Lohrasebi, Desalination 365 (2015) 176. [22] K. Sint, B. Wang, P. Kral, J. Am. Chem. Soc. 130 (2008) 16448. [23] B. Liu, R.B. Wu, J.A. Baimova, H. Wu, A.W.K. Law, S.V. Dmitriev, K. Zhou, Phys. Chem. Chem. Phys. 18 (2016) 1886. [24] Q. Tu, Q. Yang, H. Wang, S. Li, Sci. Rep. 6 (2016) 26183. [25] A.T. Nasrabadi, M. Foroutan, Desalination 277 (2011) 236. [26] D. Zhang, T. Yan, L. Shi, Z. Peng, X. Wen, J. Zhang, J. Mater. Chem. 22 (2012) 14696. [27] M. Noked, A. Soffer, D. Aurbach, J. Solid State Electrochem. 15 (2011) 1563. [28] M.E. Suss, S. Porada, X. Sun, P.M. Biesheuvel, J. Yoon, V. Presser, Energy Environ Sci. 8 (2015) 2296. [29] S. Banerjee, S. Murad, I.K. Puri, Chem. Phys. Lett. 434 (2007) 292. [30] S. Plimpton, J. Comput. Phys. 117 (1995) 1. [31] D.J. Price, C.L. Brooks, J. Chem. Phys. 121 (2004) 10096. [32] A.D. MacKerell, D. Bashford, M. Bellott, R.L. Dunbrack, J.D. Evanseck, M.J. Field, S. Fischer, J. Gao, H. Guo, S. Ha, D. Joseph-McCarthy, L. Kuchnir, K. Kuczera, F.T. K. Lau, C. Mattos, S. Michnick, T. Ngo, D.T. Nguyen, B. Prodhom, W.E. Reiher, B. Roux, M. Schlenkrich, J.C. Smith, R. Stote, J. Straub, M. Watanabe, J. Wiorkiewicz-Kuczera, D. Yin, M. Karplus, J. Chem. Phys. B 102 (1998) 3586. [33] I.S. Joung, T.E. Cheatham III, J. Chem. Phys. B 112 (2008) 9020. [34] E.M. Stuve, Chem. Phys. Lett. 519 (2012) 1. [35] G. Cicero, J.C. Grossman, E. Schwegler, F. Gygi, G. Galli, J. Am. Chem. Soc. 130 (2008) 1871. [36] J.A. Thomas, A.J.H. McGaughey, Phys. Rev. Lett. 102 (2009) 4. [37] A.A. Chialvo, P.T. Cummings, J. Phys. Chem. A 115 (2011) 5918. [38] D. Konatham, J. Yu, T.A. Ho, A. Striolo, Langmuir 29 (2013) 11884. [39] D.G. Bounds, Mol. Phys. 54 (1985) 1335. [40] H. Ohtaki, T. Radnai, Chem. Rev. 93 (1993) 1157. [41] Y. Marcus, J. Chem. Soc. Faraday Trans. 87 (1991) 2995. [42] T.A. Hilder, D. Gordon, S.H. Chung, Small 5 (2009) 2183. [43] K. Leung, S.B. Rempe, J. Comput. Theor. Nanos. 6 (2009) 1948. [44] M. Thomas, B. Corry, Philos. Trans. R. Soc. A – Math. Phys. Eng. Sci. 374 (2016) 20. [45] M. Cavallari, C. Cavazzoni, M. Ferrario, Mol. Phys. 102 (2004) 959. [46] P. Gallo, D. Corradini, M. Rovere, Phys. Chem. Chem. Phys. 13 (2011) 19814. [47] Q. Shao, J. Zhou, L.H. Lu, X.H. Lu, Y.D. Zhu, S.Y. Jiang, Nano Lett. 9 (2009) 989. [48] F. Bruni, S. Imberti, R. Mancinelli, M.A. Ricci, J. Phys. Chem. 136 (2012) 1. [49] R. Caminiti, G. Licheri, G. Piccaluga, G. Pinna, Chem. Phys. Lett. 47 (1977) 275. [50] Z.H. Duan, Z.G. Zhang, Mol. Phys. 101 (2003) 1501. [51] P. Sun, R. Ma, W. Ma, J. Wu, K. Wang, T. Sasaki, H. Zhu, NPG Asia Mater. 8 (2016) 259.