Roskilde
University
Slip length of water on graphene
Limitations of non-equilibrium molecular dynamics simulations
Kannam, S.; Todd, Billy ; Hansen, Jesper Schmidt; Daivis, Peter
Published in:
Journal of Chemical Physics
DOI:
10.1063/1.3675904
Publication date:
2012
Document Version
Publisher's PDF, also known as Version of record
Citation for published version (APA):
Kannam, S., Todd, B., Hansen, J. S., & Daivis, P. (2012). Slip length of water on graphene: Limitations of non-
equilibrium molecular dynamics simulations. Journal of Chemical Physics, 136(024705).
https://doi.org/10.1063/1.3675904
General rights
Copyright and moral rights for the publications made accessible in the public portal are retained by the authors and/or other copyright owners
and it is a condition of accessing publications that users recognise and abide by the legal requirements associated with these rights.
• Users may download and print one copy of any publication from the public portal for the purpose of private study or research.
• You may not further distribute the material or use it for any profit-making activity or commercial gain.
• You may freely distribute the URL identifying the publication in the public portal.
Take down policy
If you believe that this document breaches copyright please contact rucforsk@kb.dk providing details, and we will remove access to the work
immediately and investigate your claim.
Download date: 10. Aug. 2022
Slip length of water on graphene: Limitations of non-equilibrium
molecular dynamics simulations
Sridhar Kumar Kannam, B. D. Todd, J. S. Hansen, and Peter J. Daivis
Citation: J. Chem. Phys. 136, 024705 (2012); doi: 10.1063/1.3675904
View online: http://dx.doi.org/10.1063/1.3675904
View Table of Contents: http://jcp.aip.org/resource/1/JCPSA6/v136/i2
Published by the American Institute of Physics.
Related Articles
In situ imaging of orthoclase–aqueous solution interfaces with x-ray reflection interface microscopy
J. Appl. Phys. 110, 102211 (2011)
Interaction between water and defective silica surfaces
J. Chem. Phys. 134, 114703 (2011)
Pressure cell for investigations of solid–liquid interfaces by neutron reflectivity
Rev. Sci. Instrum. 82, 023902 (2011)
Communication: Molecular dynamics simulations of the interfacial structure of alkali metal fluoride solutions
J. Chem. Phys. 133, 061103 (2010)
Effects of nanofluids containing graphene/graphene-oxide nanosheets on critical heat flux
Appl. Phys. Lett. 97, 023103 (2010)
Additional information on J. Chem. Phys.
Journal Homepage: http://jcp.aip.org/
Journal Information: http://jcp.aip.org/about/about_the_journal
Top downloads: http://jcp.aip.org/features/most_downloaded
Information for Authors: http://jcp.aip.org/authors
Downloaded 16 Jan 2012 to 130.226.173.84. Redistribution subject to AIP license or copyright; see http://jcp.aip.org/about/rights_and_permissions
THE JOURNAL OF CHEMICAL PHYSICS 136, 024705 (2012)
Slip length of water on graphene: Limitations of non-equilibrium
molecular dynamics simulations
Sridhar Kumar Kannam,
1,a)
B. D. Todd,
1,b)
J. S. Hansen,
2,c)
and Peter J. Daivis
3,d)
1
Mathematics Discipline, Faculty of Engineering and Industrial Science and Centre for Molecular Simulation,
Swinburne University of Technology, Melbourne, Victoria 3122, Australia
2
DNRF Center ‘Glass and Time’, IMFUFA, Department of Science, Systems and Models, Roskilde University,
DK-4000, Roskilde, Denmark
3
School of Applied Sciences, RMIT University, Melbourne, Victoria 3001, Australia
(Received 2 November 2011; accepted 19 December 2011; published online 13 January 2012)
Data for the flow rate of water in carbon nanopores is widely scattered, both in experiments and
simulations. In this work, we aim at precisely quantifying the characteristic large slip length and
flow rate of water flowing in a planar graphene nanochannel. First, we quantify the slip length using
the intrinsic interfacial friction coefficient between water and graphene, which is found from equi-
librium molecular dynamics (EMD) simulations. We then calculate the flow rate and the slip length
from the streaming velocity profiles obtained using non-equilibrium molecular dynamics (NEMD)
simulations and compare with the predictions from the EMD simulations. The slip length calculated
from NEMD simulations is found to be extremely sensitive to the curvature of the velocity profile
and it possesses large statistical errors. We therefore pose the question: Can a micrometer range slip
length be reliably determined using velocity profiles obtained from NEMD simulations? Our answer
is “not practical, if not impossible” based on the analysis given as the results. In the case of high
slip systems such as water in carbon nanochannels, the EMD method results are more reliable, ac-
curate, and computationally more efficient compared to the direct NEMD method for predicting the
nanofluidic flow rate and hydrodynamic boundary condition. © 2012 American Institute of Physics.
[doi:10.1063/1.3675904]
I. INTRODUCTION
Recent developments in nanoscience and nanotechnol-
ogy have enabled the fabrication of nanofluidic devices such
as micro/nano electro-mechanical systems, nanopipettes,
nanobiosensors, nanomotors, lab-on-a-chip devices, etc.
1
As
these devices have unique attributes such as nanoliter capac-
ity, low energy dissipation, high accuracy and sensitivity, and
enhanced flow rates, understanding the physics of fluids at
the nanoscale is very important in designing, fabricating, op-
timizing, and utilizing these devices. Gravity and inertia ef-
fects which may play a dominant role at the macroscale are
negligible at the nanoscale and new phenomena due to the
high surface to volume ratio emerge.
2
Transport of momen-
tum and energy become non-local in nature
3, 4
and the molec-
ular behavior at the fluid-solid interface profoundly affects the
transport.
5, 6
Theoretically, continuum based approximations
may be invalid or difficult to formulate at the nanoscale
7
and
experimental difficulties are also not yet fully solved.
8
Com-
puter simulations of the type used in the present study have
therefore become a viable tool for studying nanoscale phe-
nomena.
Starting from the last decade, water confined in carbon
nanostructures has received significant attention due to the
importance of water and the unique properties of carbon.
9–43
a)
Electronic mail: urssrisri@gmail.com.
b)
Electronic mail: btodd@swin.edu.au.
c)
Electronic mail: jschmidt@ruc.dk.
d)
Electronic mail: peter.daivis@rmit.edu.au.
A number of studies have been aimed at quantifying the slip
length and flow rate of water in carbon nanotubes (CNTs)
and flat graphene nanochannels, both experimentally
12–19
and
in simulations.
21–38
However, the data are widely scattered:
even over orders of magnitude and no consensus has been
reached. Water has been shown to have almost a plug like
velocity profile, resulting in 1 to 5 orders of magnitude flow
enhancement compared to classical Navier-Stokes prediction
assuming a no-slip boundary condition (BC). Qualitatively,
some studies have found a very high slip length in very small
diameter CNTs and as the tube diameter increases, the slip
length approaches a constant value which is equal to the slip
length on a flat graphene surface.
29, 32, 33, 43
Contrary to this,
some researchers have found the opposite behavior and at-
tribute this to the increase in surface friction as the tube diam-
eter is decreased.
23–25
In a very recent field effect transistor
experimental study, Qin et al.
18
found non-monotonic behav-
ior of slip, which has also been reported by Sokhan et al.
44
for
slip of methane in CNTs. Hence, even qualitatively we have
contradictory slip behavior in CNTs.
In experimental studies, Maali et al.
17
found a slip length
of 8 ± 2 nm for water on a graphite surface. Qin et al.
18
found a decreasing slip length as the CNT diameter is in-
creased from 0.81 to 1.59 nm with a non-monotonic behavior
at 1.08 nm diameter tube. As the tube diameter increases to
1.59 nm, its slip length converges to 10 nm, which can be ap-
proximated to the slip length of water on planar graphene. Zhu
and Granick
45
and Tretheway and Meinhart
46
found microm-
eter range slip lengths for water on hydrophobic surfaces. In
0021-9606/2012/136(2)/024705/9/$30.00 © 2012 American Institute of Physics136, 024705-1
Downloaded 16 Jan 2012 to 130.226.173.84. Redistribution subject to AIP license or copyright; see http://jcp.aip.org/about/rights_and_permissions
024705-2 Kannam
et al.
J. Chem. Phys. 136, 024705 (2012)
simulations, Thomas and McGaughey
29
found 30 nm, Falk
et al.
32
found 80 nm, and Kotsalis
34
found 67 ± 45 nm slip
lengths for water on graphene surfaces. Babu and Sathian
33
found monotonically decreasing slip length as the tube di-
ameter is increased from 0.81 to 5.42 nm and in the higher
diameter tubes the slip length converges to just about 1 nm,
which can be taken as the slip length on a flat graphene sur-
face. Some of this variation can be attributed to the difference
in water models and carbon models used in simulations. The-
oretically, Myers
43
predicted 39 nm slip length in the limit of
high diameter CNTs. As explained by Kotsalis
34
the large un-
certainties (67 ± 45 nm) are due to the large variation in the
velocity gradient of water at the graphene surface. Thomas
and McGaughey
29
calculated the viscosity of confined wa-
ter in equilibrium molecular dynamics (EMD) simulations by
using the Einstein self-diffusion coefficient and they used this
viscosity to constrain the fit to the non-equilibrium molecular
dynamics (NEMD) streaming velocity profiles. They found
that unconstrained fits result in almost 100% variation in the
predicted slip length and the viscosity of water.
Experimentally, fabricating a perfect defect free individ-
ual carbon nanotube is extremely difficult. Moreover measur-
ing the tube diameter, controlling the pressure difference to
drive the fluid, and finally performing the nanoliter volume
experiment is a cumbersome procedure, which could be one
reason for the scattered data in experiments.
12–19
In simula-
tions, the commonly used NEMD methods also have their
limitations. At room temperature water has an average ther-
mal velocity of 340 m/s and the fluid velocities in experiments
are on the order of 0.01 m/s. NEMD simulations can only
simulate the systems for a few nanoseconds with a time step
of ∼1 fs. Due to this computational limitation, NEMD sim-
ulations are done with very high pressure gradients (or shear
rates in Couette flow) to obtain a mean fluid velocity compa-
rable to the thermal velocity and thus statistically significant
results. At these high fields nonlinear effects may begin to
emerge and the slip length diverges.
47
Therefore, the extrap-
olation of NEMD results to experimental fields is not reliable
and is likely to lead to deviation from the flux determined un-
der experimental conditions, where the hydrodynamic proper-
ties obey linear relations.
41
Moreover, to do the extrapolation
these NEMD simulations have to be performed for a range of
pressure gradients or shear rates. In order to reduce the com-
putational time and thus avoid the complex carbon models,
a number of studies freeze the carbon atoms to their lattice
sites and thermostat the water directly to maintain the desired
temperature. It has been shown that the molecular momentum
transfer at the fluid-solid interface plays a key role in nanoflu-
idic behavior as the fluid transport properties are dominated
by the interface.
5, 44, 48
Indeed, thermostatting the fluid instead
of the walls can lead to unphysical distortions of the fluid’s
stress response and can artificially enhance or reduce the slip,
depending on the type of wall surface.
44, 48
Sokhan et al.
44, 49
observed an increase of ∼20% in slip for methane in flexi-
ble graphene/CNTs compared to rigid counterparts. Refer to
our previous work for a broader discussion.
47
It is however
worth pointing out here that a systematic study of the ef-
fects of wall roughness, frozen vs vibrating walls (and cor-
respondingly, thermostatting the fluid vs thermostatting the
walls), and the various ab initio slip models is yet to be ac-
complished and would represent a very useful study in the
field. For example, the recent publication by Groombridge
et al.
37
uses the equilibrium model of Sokhan and Quirke
50
to study the slip of water confined by amorphous surfaces of
carbon and polydimethylsiloxane and finds reasonable agree-
ment between NEMD and EMD results (see their Table 1).
For high slip systems such as water against a hydropho-
bic surface, the problem becomes even more complicated. As
mentioned above, and as we will show in Sec. III,avery
small change in the measured NEMD velocity profile can
result in a very large deviation in the slip length. This sug-
gests the need for developing new theoretical methods for the
prediction of such a highly sensitive phenomenon. Recently,
Hansen et al.
51
proposed a method of calculating the intrin-
sic friction between fluid and solid using EMD simulations
based on a statistical mechanical approach. The method has
successfully predicted the slip length for atomic fluids and
methane confined between Lennard-Jones solid walls
51
and
inside graphene nanochannels.
47
Our preliminary results of
methane flow in CNTs show that the minimum slip length
for any diameter CNT is greater than or at least equal to the
slip length on a flat graphene surface.
52
Understanding the
behavior of water on a flat graphene surface and documenting
its slip length is crucial before we completely understand the
behavior in CNTs which are a cylindrical form of graphene,
where the curvature of the tube affects the friction. This is
the motivation of the present study. In this work using EMD
simulations (with no pressure gradient or shear rate), thus by-
passing limitations of NEMD methods, we calculate the inter-
facial friction coefficient between water and graphene. Using
this friction coefficient we determine the slip length and flow
rate. We also perform NEMD simulations of Poiseuille and
Couette flows for a wide range of external fields and shear
rates to understand the slip phenomena further and compare
these results to our EMD method predictions.
II. SIMULATION DETAILS
We use the recently parameterized flexible simple point
charge (SPC/Fw) model for water, which reproduces the dy-
namical properties of water close to experimental values.
53–55
Graphene is modeled using the second generation reactive
empirical bond order Tersoff-Brenner potential, which is
widely used for carbon allotropes.
56
Electrostatic interac-
tions between water molecules are modeled using the Wolf
method,
29, 30, 55, 57, 58
which enables us to simulate for longer
times and 20 simulations at each state point. The interaction
between water molecules and carbon atoms of the graphene
is modeled using the Lennard-Jones potential with parame-
ters of Werder et al.
42
We have used two layers of graphene
at two walls to produce a stable system and better heat
conduction between fluid and solid. A weak Lennard-Jones
potential is applied between carbon atoms belonging to the
different graphene layers to hold them together at 0.34 nm
distance. All Lennard-Jones interactions are truncated at a dis-
tance of 1 nm. The channel width, i.e., the distance between
the two innermost graphene layers is set to 3.9 nm (roughly 12
molecular diameters) in the y direction and periodic boundary
Downloaded 16 Jan 2012 to 130.226.173.84. Redistribution subject to AIP license or copyright; see http://jcp.aip.org/about/rights_and_permissions
024705-3 Water-graphene J. Chem. Phys. 136, 024705 (2012)
conditions are applied along the x and z directions. The van
der Waals size of the carbon atoms (0.34 nm) is subtracted
from 3.9 nm to define the available channel width. Room tem-
perature (300 K) is maintained by applying the Nosé-Hoover
thermostat to carbon atoms, so that the heat produced by
the slip and viscous friction in the water is conducted away
through the graphene as is done in a real experiment. The wa-
ter density in the channel is kept at 1000 kg/m
3
. Simulation
time spans are from 5 to 10 ns with a time step of ∼1fsusing
the leapfrog integration algorithm.
III. RESULTS
A. EMD simulations
We refer to Hansen et al.
51
for details on the method of
calculating the friction coefficient. In brief, we choose a water
slab of width one molecular diameter close to the graphene
wall.
47, 59
After equilibration, we evaluate the x direction wall-
slab shearing force via
F
x
(t) =
i∈slab
j∈wall
F
ij,x
(t)(1)
and the centre of mass (CM) velocity of the slab via
u
slab
(t) =
1
m
i∈slab
m
i
v
i,x
(t), (2)
where F
ij, x
is the force on a slab water molecule i due to the
carbon atom j at time t (the carbon-oxygen van der Waals
force), v
i, x
is the x-component of the velocity of the slab water
molecule i and m is the total mass of the slab m =
i ∈ slab
m
i
.
Using these two quantities, we evaluate the slab velocity-force
cross correlation function C
uF
x
(t) and slab velocity autocor-
relation functions C
uu
(t)
C
uF
x
(t)=u
slab
(0)F
x
(t) and C
uu
(t) =u
slab
(0)u
slab
(t).
(3)
The Laplace transforms,
C
uF
x
(s) and
C
uu
(s), are then related
via the constitutive equation
51
C
uF
x
(s) =−
ζ (s)
C
uu
(s) . (4)
The zero frequency coefficient, ζ
0
, is found via the fit be-
tween the above two Laplace transformed correlation func-
tions assuming that
ζ (s) is a Maxwellian memory function.
Finally, the friction coefficient, ξ
0
, is calculated by dividing
ζ
0
by the surface area of the graphene A. We again refer to
Hansen et al.
51
for full details of the steps involved. Our re-
sult for the friction coefficient between water and graphene is
(1.25 ± 0.10) × 10
4
kg m
−2
s
−1
.
The Navier
60
slip length is defined as
L
s
=
η
0
ξ
0
. (5)
Using the bulk water shear viscosity η
0
= (7.5 ± 0.5)
× 10
−4
kg m
−1
s
−1
, the slip length of water on a planar
graphene surface is thus estimated as 60 ± 6nm.
B. NEMD simulations
We perform both Poiseuille and Couette flow NEMD
simulations for a wide range of external fields and shear rates
starting from the lowest possible values. Poiseuille flow is
generated by applying a constant external field to all the atoms
of the water molecules and Couette flow is generated by mov-
ing the upper graphene wall with a constant velocity while
keeping the lower graphene wall fixed. We fit the Poiseuille
flow velocity profiles to a quadratic equation u
x
(y) = ay
2
+ b and Couette flow velocity profiles to a linear equation
u
x
(y) = ay + b (see Fig. 1). Using these fits one can deter-
mine the slip velocity and the fluid velocity gradient at the
wall. From these two quantities the NEMD slip length can be
found via
u
s
= L
s
∂u
x
∂y
y=y
w
, (6)
see Fig. 1. The slip length and friction coefficient are both
intrinsic properties of the fluid-solid interface and are inde-
pendent of the flow type.
1. Couette flow
To our knowledge Couette flow type simulations have not
previously been carried out on the water-graphene system in
order to determine the slip length. We move the upper wall
with a constant velocity in the range 5 to 1000 m/s. Within
this range, linearity is expected to hold over a small range of
low fluid velocities, above which the slip length is no longer
a constant. Even the lowest applied shear rate (wall velocity
divided by channel width) is 3 to 4 orders of magnitude higher
than the highest experimental shear rate.
61
When we have a
small slip system, we can use low wall velocities as the wall
FIG. 1. Definitions of slip length and slip velocity for Poiseuille and Couette flows.
Downloaded 16 Jan 2012 to 130.226.173.84. Redistribution subject to AIP license or copyright; see http://jcp.aip.org/about/rights_and_permissions
