University of Southern California Dissertations and Theses

Dynamics of water in nanotubes: liquid below freezing point and ice-like near boiling point

(USC Thesis Other)

Dynamics of water in nanotubes: liquid below freezing point and ice-like near boiling point

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content Dynamics of Water in Nanotubes: Liquid Below Freezing Point and Ice-like Near
Boiling Point
by
Jos e Cobe~ na-Reyes
A Dissertation Presented to the
FACULTY OF THE GRADUATE SCHOOL
UNIVERSITY OF SOUTHERN CALIFORNIA
In Partial Fulllment of the
Requirements for the Degree
DOCTOR OF PHILOSOPHY
(CHEMICAL ENGINEERING)
May 2021
Copyright 2021 Jos e Cobe~ na-Reyes
Acknowledgements
Finishing something as complex as a Ph.D. program requires help from dierent sources. I feel
fortunate that I was able to nd many extraordinary people along the way in this journey. Without
these people it would have been dicult to cross the nish line. In the following paragraphs I
acknowledge them.
First, I would like to thank to my advisor, Professor Muhammad Sahimi. His patience, exper-
tise and guidance has helped me to become an independent researcher. Professor Sahimi always
set high standards and led by example, being ambitious, working hard and never being afraid to
tackle new challenges in our scientic elds. The ndings presented in this thesis compile new
advances regarding water under connement in nanotubes of several types that we hope will be
helpful to others. I am grateful to have been a part of his research team.
I would like to thank my defense committee, Professors Rajiv Kalia and Aiichiro Nakano. My
rst paper was published thanks to Professor Sahimi and Professor Kalia who gave me important
feedback since my qualication exam. I consider Professor Nakano a co-advisor. I took a few
classes with him and he always pushed me to do my best, without fear of failure and to be eager
to learn more, to be bold and to be ambitious.
I have been lucky to be a part of the Mork Family Department. Here, I met Andy Chen,
who was rst the Director of Student Aairs, and then became the Business Manager of the
department. Now, having him as the Director of Doctoral Programs of the Viterbi School of
Engineering has been extremely helpful. Andy has always showed empathy and patience, caring
about the students and nding solutions to any issues we faced. Similarly, I want to thank to
ii
Maria Albertina (aka Tina) Silva who welcomed me into the instructional laboratories. Thanks
to Tina and Professor Sahimi I was able to land TA positions several times. It was always nice to
run into her inside the building and casually chat for a few minutes. I am grateful for their help.
My dearest SHPE-USC friends Joel, Sophia, Spencer and Emily. Who made my journey a
little more pleasant and enjoyable. I consider them an important part of my support system
during the Ph.D. I can't imagine how my life would have been by now without their help.
I would also like to thank to Diana Ramos. Her love and patience were signicant to accomplish
this. Diana will always have a special place in my mind and heart.
I want to thank my friends and oce-mates Nariman, Christine and Zumra. We shared many
stories and little adventures that will always be in my memory. Thank for you making my stay a
little better.
I also want to thank to my high school friends, or as they call themselves Miserables. Their
[virtual] company was enjoyable, especially during these pandemic times. Similarly, I would like
to thank to my college friends Paula and David, who always encouraged me to do my best.
As a Ph.D. student I realized that my circle of friends constantly changed. People like Sara,
Laura, Debora, Lorena and Sebastian were present at dierent stages and helped me have a more
pleasurable time in Los Angeles. I am grateful for the time we spent together.
There were many people that helped me in dierent ways. I want to thank to Professor Andrea
Hodge, Chair of the Mork Family Department and Dr. Brandi Jones, Vice Dean for Diversity
and Strategic Initiatives for helpful talks and advice. I also want to thank to the Center of
Engineering Diversity (CED) and their director Tracy Navarro for creating safe spaces oriented
to minorities and underrepresented groups in engineering. Additionally, I would like to thank
Professors Robert Young, Jincay Chang, Katherine Shing and Thomas Marlin for allowing me be
part of their teaching team. I would also like to thank to Broderick, Parissa and Sergio. Your
help can't be measured. You helped me to become a better human being.
iii
There were people that were not directly involved during the time I was in the Ph.D. program
but who helped me pave the path and whose impact in my life was instrumental in my success.
My rst advisor, who left us too soon, Professor Carlos Mu~ noz Cajiao, who was my undergraduate
thesis advisor and a good friend. He was the rst Professor to show me how to model and simulate
processes. I am eternally grateful for everything he did for me. I also want to thank to one of
the smartest Professors I have ever known, Professor William Loyola, who was my M.B.A. thesis
advisor. I want to thank him for his talks and wise advice.
I also want to thank to SENESCYT from Ecuador and the Petroleum Research Fund for their
nancial support.
Above all, I would like to thank to my family. My mother, Isabel Reyes, and my father, Jos e
Cobe~ na. Their strength and patience have always motivated me to continue despite the obstacles
I found in the way. My siblings Andr es and Isabel and my brother-in-law Danilo, I am grateful
for your strength, support and patience. Your help was crucial to my success.
iv
Table of Contents
Acknowledgements ii
List of Figures vii
Abstract ix
Chapter 1: Introduction 1
1.1 Background and motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.2 General structure of nanotubes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.3 Structure, properties and uses of nanotubes . . . . . . . . . . . . . . . . . . . . . . 4
1.3.1 Carbon Nanotubes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.3.2 Silicon-Carbide Nanotubes . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.4 Water conned in nanotubes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
1.5 Computational Techniques and Models . . . . . . . . . . . . . . . . . . . . . . . . . 9
1.5.1 Molecular dynamics simulation . . . . . . . . . . . . . . . . . . . . . . . . . 9
1.5.1.1 Molecular interactions and Hamiltonian mechanics . . . . . . . . . 10
1.5.1.2 Velocity Verlet algorithm . . . . . . . . . . . . . . . . . . . . . . . 11
1.5.2 Coulombic interactions and long-range summation . . . . . . . . . . . . . . 12
1.5.3 Force Fields . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13
1.5.3.1 TIP4P water models . . . . . . . . . . . . . . . . . . . . . . . . . . 13
1.5.3.2 Bond Order Potentials: Terso and REBO . . . . . . . . . . . . . 15
1.5.3.3 Reactive Force Field: RexPoN . . . . . . . . . . . . . . . . . . . . 17
1.5.4 Monte Carlo simulation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
Chapter 2: Complex Behavior of Ordered and Icelike Water in Carbon Nanotubes
Near its Bulk Boiling Point 19
2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
2.2 Computational Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21
2.2.1 Molecular Models and MD Simulation Protocol. . . . . . . . . . . . . . . . 21
2.2.2 Kirkwood g-Factor. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
2.2.3 The ten Wolde Parameter. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
2.2.4 Mean-Square Displacement . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
2.2.5 Radial Distribution Function. . . . . . . . . . . . . . . . . . . . . . . . . . . 25
2.2.6 Mean Connectivity. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
2.2.7 Cage Correlation Function. . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
2.2.8 Intermediate Scattering Function. . . . . . . . . . . . . . . . . . . . . . . . 27
2.3 Results and discussions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
2.3.1 The Kirkwood g-Factor. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
2.3.2 The ten Wolde Parameter. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
2.3.3 Dependence of the Potential Energy on the Radial Distance. . . . . . . . . 31
2.3.4 The Cage Correlation Function. . . . . . . . . . . . . . . . . . . . . . . . . . 32
v
2.3.5 Intermediate Scattering Function. . . . . . . . . . . . . . . . . . . . . . . . 33
2.3.6 Mean-Square Displacements. . . . . . . . . . . . . . . . . . . . . . . . . . . 34
2.3.7 Radial Distribution Function and Mean Connectivity of the Atoms. . . . . 35
2.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37
Chapter 3: Rheology of water in small nanotubes 38
3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
3.2 Molecular Models and Computational Methods . . . . . . . . . . . . . . . . . . . . 42
3.2.1 The Nanotubes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42
3.2.2 Computational protocol for bulk water . . . . . . . . . . . . . . . . . . . . . 43
3.2.3 Computational protocol for conned water . . . . . . . . . . . . . . . . . . 45
3.3 Results and discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
3.3.1 The cage correlation function and mean-squared displacements . . . . . . . 50
3.3.2 Radial distribution function and density
uctuations of liquid water . . . . 51
3.3.3 Intermediate scattering function and the relaxation time . . . . . . . . . . . 54
3.3.4 Rheology . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
3.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
Chapter 4: Universal Intrinsic Dynamics and Freezing of Water in Small Nan-
otubes 61
4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
4.2 Molecular models and details of the computations . . . . . . . . . . . . . . . . . . 65
4.2.1 Computational models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
4.2.2 Dipole Orientation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67
4.2.3 The Axial and Radial Distribution Functions . . . . . . . . . . . . . . . . . 67
4.2.4 The Diusion Coecient . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68
4.2.5 The Cage Correlation Function . . . . . . . . . . . . . . . . . . . . . . . . . 68
4.2.6 Hydrogen Bond Correlation Function . . . . . . . . . . . . . . . . . . . . . 70
4.3 Results and discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
4.3.1 Dipole orientation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
4.3.2 Ice types . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72
4.3.3 The Radial and Axial Distribution Functions . . . . . . . . . . . . . . . . . 73
4.3.4 The Diusion Coecient . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75
4.3.5 The Cage Correlation Function . . . . . . . . . . . . . . . . . . . . . . . . . 77
4.3.6 Strength of Hydrogen-Bond Networks . . . . . . . . . . . . . . . . . . . . . 77
4.3.7 Eect of Partial Charges. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79
4.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80
Chapter 5: Chemical potential and phase change of water in CNTs using RexPoN
potential 81
5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81
5.2 Computational Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83
5.2.1 Description of the RexPoN potential . . . . . . . . . . . . . . . . . . . . . . 83
5.2.2 Chemical potential of the RexPoN water model . . . . . . . . . . . . . . . . 84
5.2.3 Filling up the CNTs with water . . . . . . . . . . . . . . . . . . . . . . . . . 86
5.2.4 Molecular dynamics of water in CNTs using RexPoN . . . . . . . . . . . . . 87
5.3 Results and discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87
5.4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90
Bibliography 91
vi
List of Figures
1.1 Representation of the vectors a
1
, a
2
, and C
h
in a nanotube. In this example the
lattice translational indices are n 8 and m 4. . . . . . . . . . . . . . . . . . . . 4
2.1 The Kirkwood g-factor calculated using periodic boundary conditions in all direc-
tions and the PPPM method for computing the Coulombic interactions. . . . . . . 29
2.2 The Kirkwood g-factor at several temperatures. . . . . . . . . . . . . . . . . . . . . 29
2.3 Comparison of the ten Wolde parameter for the bulk and conned water atT 363
K. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
2.4 Snapshot of the distribution of water inside the carbon nanotube, with oxygen
(red), hydrogen (white), and hydrogen bonds (light blue lines). . . . . . . . . . . . 31
2.5 Dependence of the potential energy Uprq on the radial distance r from the CNT's
center. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
2.6 The cage correlation function Cptq. (a) The overall function at several tempera-
tures. (b) Comparison of Cptq for the water molecules near the tube's center and
in the outer layer. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
2.7 The intermediate scattering function at 343 K. . . . . . . . . . . . . . . . . . . . . 34
2.8 Axial mean-square displacementsxR
2
ptqy of water (oxygen) molecules. . . . . . . . 35
2.9 The radial distribution function gprq at 363 K. . . . . . . . . . . . . . . . . . . . . 36
3.1 The cage correlation function Cptq. The results at 80 K indicate a frozen state,
whereas those at other temperatures are representative of a liquid state. . . . . . . 50
3.2 Axial mean-squared displacementsx
2
zptqy of water. The results at 80 K indicate
a frozen state. The inset compares the diusivities in connement with the bulk
values. The bulk values were taken from Weiss et al. (2011) . . . . . . . . . . . . . 51
3.3 The radial distribution functiongprq, which provide evidence for a (a) liquid state,
and the existence of (b) low- and (c) high-density water regions in the (12,0) nanotube 53
3.4 Histogram of density
uctuations of water in the (12,0) nanotube . . . . . . . . . . 54
vii
3.5 The intermediate scattering function F
s
p;tq for (a) bulk water, and (b) water
conned in the nanotube. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
3.6 Viscosity of bulk water as a function of the shear rate, computed by the NEMD
simulation. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56
3.7 Zero-shear rate viscosity of water in the nanotube, computed by the Green-Kubo
equation. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56
3.8 Temperature- and shear rate-dependence of the shear stress in the nanotube, com-
puted by the NEMD simulation. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
3.9 Temperature- and shear-rate dependence of the eective viscosity of water, com-
puted by the NEMD simulation. They indicate non-Newtonian rheology. . . . . . . 58
4.1 Histogram of the dipole angle of water conned in CNTs. The highest peaks
indicate the most frequent dipole angle, which is around 30

. . . . . . . . . . . . . 71
4.2 Histogram of the dipole angle of water conned in SiCNTs. The sharp drops in
the peaks indicate the possible melting points. The most frequent dipole angle is
close to 40

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72
4.3 The structure of the various types of ice in the SiCNTs. . . . . . . . . . . . . . . . 73
4.4 The axial distribution function of water in the SiCNTs. Similar to Fig. 4.1, a
sudden decrease in the peaks indicate a rst-order transition. . . . . . . . . . . . . 74
4.5 The radial distribution function of water in the SiCNTs. Similar to Figs. 4.2 and
4.4, a sudden decrease in the peaks indicate a rst-order transition. . . . . . . . . . 75
4.6 Temperature-dependence of the axial diusion coecientD
z
, using the VACF (left)
and the MSDs (right). Arrows at lowerT indicate the possible melting point, while
those at higher T indicate the fragile-to-strong transition. . . . . . . . . . . . . . . 76
4.7 Cage correlation function Cptq of the oxygen atoms. The decay of Cptq occurs at
temperatures that are consistent with those in Figs. 4.2, 4.4 and 4.5. . . . . . . . . 77
4.8 Hydrogen-bond correlation function C
HB
of water in SiCNTs (circles) and CNTs
(squares). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78
4.9 Eect of (articially) placed partial charges on the hydrogen bond correlation func-
tion C
HB
of water in the CNTs at T 273 K. The rst two upper curves are the
same as those in Fig. 4.8. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79
5.1 Excess chemical potential of water . . . . . . . . . . . . . . . . . . . . . . . . . . . 88
5.2 Ice phases observed in the three CNTs. . . . . . . . . . . . . . . . . . . . . . . . . . 88
5.3 Potential energy of water in the three CNTs . . . . . . . . . . . . . . . . . . . . . . 89
5.4 Radial distribution function gprq of water below and above melting point. . . . . . 89
viii
Abstract
Water under connement is found everywhere: Inside fractured rock, within blood vessels and
inside living tissue and cells. Nonetheless, the behavior of bulk and conned water still puzzles
scientists around the world. There are many questions that have not been completely answered
yet.
The rise of the nanotechnology brought new materials with many potential applications. Car-
bon and Silicon Carbide Nanotubes (CNTs and SiCNTs) showed remarkable properties such as
very high Young's modulus, extremely high thermal conductivity, capacity to be tuned and doped,
and the potential to behave as metallic, non-metallic and semiconductors materials. Since the
discovery of nanotubes two decades ago, researchers started the dicult task to study these ma-
terials and how the properties of water departs from the bulk when it is conned inside these
nanotubes.
Properties of water such as the melting and boiling points, viscosity, or even ice phases can
suer drastic changes when conned inside nanotubes. It is the purpose of this Thesis to identify
and understand such changes. To achieve this goal we employed several computational techniques
and models. Molecular Dynamics (MD) and Monte Carlo (MC) simulations were intensively used
in this work. Similarly, accurate force elds such as the REBO potential and the RexPoN potential
were also employed.
Our methods and results are organized in this Thesis as follows: In Chapter 1, a brief explana-
tion of the computational methods used in this Thesis is presented. Chapter 2 shows the results
of our rst work where we used water conned in an isolated CNT and showed that water can
ix
exhibit an ice-like behavior even when it reaches temperatures close to the bulk boiling point. In
Chapter 3 we present how the connement can change the rheology of water. We show that a
non-Newtonian behavior appears in water when it is conned in a SiCNT. The non-Newtonian
behavior appears \early" compared with bulk water. In Chapter 4 we seek to understand whether
the changes in melting point when water is conned, are universal. We also studied the eect of
the partial charges in the change of the melting point. Finally, in Chapter 5 we used the RexPoN
potential to compute the excess chemical potential to use it in Monte Carlo simulations to create
structures to be used in MD simulations. We found that there is a change in their melting point,
but not as dramatic as initially though. The purpose of Chapter 5 is to construct a phase diagram
using the RexPoN force eld.
We hope that the ndings presented in this work will be useful to the scientic community
and that nano
uidic systems can be developed based on our results.
x
Chapter 1
Introduction
1.1 Background and motivation
Water in conned media has a unique behavior. Its dynamic and static properties widely dif-
fer from those of bulk water. The signicance of this behavior has motivated enormous research
eorts, as understanding water properties in connement is relevant to a wide variety of prob-
lems. Perhaps most importantly, conned water is commonly encountered in nature, such as, for
example, inside cells, in blood vessels, or in fractured rock. The study of such conned water is
relevant to many applications, starting with the fabrication of desalination devices, drug-delivery
and biomedical devices, gas storage, and fuel cells, among others.
Nonetheless, since the discovery of C60 structures by Kroto et al. (1985) and later, in 1991,
the discovery of carbon nanotubes (CNTs) by Iijima (1991), CNTs and other similar materials,
such as silicon-carbide nanotubes (SiCNTs) and boron-nitride Nanotubes (BNNTs) have received
an extraordinary amount of attention in the last several years, not only because of the unusual
properties that these materials exhibit, but also due to the changes that water undergoes when it
is conned inside nanotubes. Several studies, both experimental (Bertrand et al., 2013b; Cupane
et al., 2014; Diallo et al., 2015) and computational (Chen et al., 1997; Brovchenko et al., 2003;
Debenedetti, 2003; Alba-Simionesco et al., 2006; Banys et al., 2006; Bertrand et al., 2013a; Zhao
1
et al., 2010) have been performed by scientists in order to understand, for example, how phase
transitions change (Takaiwa et al., 2008) in nanotubes; how the mechanical (Liu and Wang, 2005;
Zhang et al., 2011; Kohler and da Silva, 2016) and transport properties are altered (Liu et al.,
2005) and, in general, how water behavior in conned media departs from the bulk behavior.
More specically, recent experiments performed by Agrawal et al. (2017a) and Chiashi et al.
(2019), as well as several molecular dynamics simulations studies (Koga et al., 2001; Takaiwa
et al., 2008; Kyakuno et al., 2011) on water in CNTs have revealed that ice-like water structure
can be found at temperatures that are dierent from the bulk melting point. Similarly, Khademi
and Sahimi (2016) showed that water conned in a SiCNT actually does not freeze, even at
temperatures as low as 100 K. Such discoveries raise the question of whether water can freeze
inside any type of nanotube. They also give rise to questions regarding the transport properties
of water when it does not freeze and remains as supercooled water.
This dissertation aims to address these problems, including the question of ice-like behavior of
water in CNTs even at the bulk boiling temperature. In addition, if supercooled water in certain
nanotubes does not freeze, then the question arises as to whether such state of water still behaves
as a Newtonian
uid, similar to bulk water, or whether it has a more complex non-Newtonian
behavior.
In this chapter we rst present a general overview of nanotubes used in this research. Then,
the methods employed are described in the rest of the chapter.
Chapters 2, 3 and 4 have already been published. The reader can nd the published versions
here: Chapter 2: Cobe~ na Reyes et al. (2018), Chapter 3: Cobe~ na-Reyes and Sahimi (2020a) and
Chapter 4: Cobe~ na-Reyes and Sahimi (2020b).
2
1.2 General structure of nanotubes
Nanotubes are cylindrical and hollow structures that contain one or more type of atoms ar-
ranged in the body of the tube. Their structure can be described by two integers, pn;mq, which
are the lattice translational indices. These values are used to form what is called the chiral vector
C
h
. The magnitude of C
h
also depends on the values of a
1
and a
2
, which are the unit vectors
representing the unit cell:
C
h
na
1
ma
2
(1.1)
The unit vectors a
1
and a
2
are dened in terms of the bond length a
bl
:
a
1

C
3
2
a
bl
;
?
3
2
a
bl
G
(1.2)
a
2

C
3
2
a
bl
;
?
3
2
a
bl
G
(1.3)
From C
h
the diameter of the tube d
t
is calculated using the values ofpm;nq and a
bl
:
d
t

?
3a
bl
a
pn
2
nmm
2
q

(1.4)
The chiral angle is dened as the angle between the vectors C
h
and a
1
:
cos
2nm
2
?
n
2
nmm
2
(1.5)
When 30

ornm, the nanotube is calledarmchair. If 0 orm 0, thepn; 0q nanotube
is called zigzag. If 0 30 or nm and n 0, the tube is called chiral. Figure 1.1 shows a
nanotube with n 8 and m 4.
3
Figure 1.1: Representation of the vectors a
1
, a
2
, and C
h
in a nanotube. In this example the
lattice translational indices are n 8 and m 4.
Depending on the values of pn;mq, there are nanotubes that are more stable than others.
This situation depends on the internal interactions and forces such as the partial charges and
electrostatic forces. In the case of CNTs, the partial charges are considered to be 0. Additionally,
a nanotube is calledsinglewalled (SW) if it is composed of only one single cylinder; doublewalled
(DW) if it is made of two concentric cylinder and, more generally, one may speak of multiwalled
(MW) nanotubes, when the number of concentric tubes is grater than two. In this work, however,
only singlewall nanotubes have been studied
1.3 Structure, properties and uses of nanotubes
The two types of nanotubes (CNTs and SiCNTs) that we used in this work share similar
properties. In this section their structure, applications and characteristics are explained.
4
1.3.1 Carbon Nanotubes
CNTs are a form of carbon with cylindrical shape, and belong to the family of fullerene
structures. The diameter of the tubes can be as small as 4

A(Qin et al., 2000), up to a few
nanometers. Since they are made of only carbon, the C atoms are expected to have 0 or negligible
partial charges. The length of the C-C bond is 1:42

A, according to the calculations by density
functional theory (DFT) (Budyka et al., 2005).
CNTs exhibit outstanding mechanical, thermal and electrical properties. They have very high
Young's modulus, high stiens, and high axial strength (Treacy et al., 1996; Coleman et al., 2006).
The Young's modulus reported is higher than 1 TPa (Lau et al., 2006), higher than diamond's
Young's modulus. Similarly, it has been shown that the thermal conductivity of CNTs is around
6600 W/mK (Berber et al., 2000), which means that they are excellent thermal conductors. Their
electronic properties, on the other and, depends on the chirality. Zigzag CNTs exhibit metallic
properties and can be very good electric conductors (Odom et al., 1998), whereas chiral and
armchair CNTs can behave as metallic or semiconducting materials, depending on the pn;mq
values.
CNTs have found a place in composite materials. Due to their high stiness and Young's
modulus, CNTs can provide extra strength and enhance the performance of other materials. For
example, CNTs composites made with polymeric matrices have been proposed for fabricating
aircraft and aerospace parts (Gohardani et al., 2014). Additionally, metallic and metalloid com-
posites for pure metals and their alloys, such as Sn-Ag-C-CNT, Al-CNT, Cu-CNT and Si-CNT,
have been tested as solder materials, in electronic packaging, and as parts in Li-ions batteries,
respectively (Agarwal et al., 2016).
CNTs have been also shown to be excellent adsorbents, since it was found that they can
be used to remove heavy metals, such as Zn, Cd (Abbas et al., 2016) and Pb (Saleh, 2016) from
water. Applications related with energy storage and distribution have been also discussed (Gu and
5
Yushin, 2014; Lee et al., 2014; Wen et al., 2016), due to CNTs' electrical properties. Likewise,
biomedical applications, such as cancer detection, implants manufacturing, and drug and gene
delivery have been the subject of recent studies (Eatemadi et al., 2014; Martincic and Tobias,
2015; Karimi et al., 2015). Finally, CNTs are also good candidates for storage of gas natural
(Mahdizadeh and Goharshadi, 2013) and hydrogen (Jin et al., 2015); they could also be used
to detect such gases as NH
3
and CO
2
(Abdelhalim et al., 2015), also as part of sports goods
(Kingston et al., 2014) and in cement composites (D'Alessandro et al., 2016).
1.3.2 Silicon-Carbide Nanotubes
As their name suggests, SiCNTs are made of of Si and C atoms. Theoretically, many combi-
nations of Si and C atoms could exist. It has been shown, however, that the SiCNTs that have
the same number of atoms of C and Si with a zigzag conguration are the most stable in terms
of their internal energy (Menon et al., 2004). A particular characteristic of SiCNTs is their bond
length, which is estimated to be 1.80

A, longer than the C-C bond in CNTs. Moreover, since
both atoms are dierent in SiCNTs, partial charges are expected to arise. It was shown, based on
DFT calculations, that the partial charges for SiCNTs with the same number of atoms of Si and
C are +0.6e and -0.6e, respectively (Menon et al., 2004). On the other hand, while CNTs can be
either metallic or semmiconducting, SiCNTs are always semiconductors.
SiCNTs have similar applications than CNTs. In addition, they have similar values of Young's
modulus, around 0.6 TPa. (Petrushenko and Petrushenko, 2015). Nonetheless, compared with
CNTs, SiCNTs have a higher resistant to corrosive environment (Khademi and Sahimi, 2011).
Moreover, SiCNTs have been employed in the fabrication of membranes for gas separation (Elyassi
et al., 2007), and can be used in environments at very high temperatures (Khademi and Sahimi,
2011). Additionally, molecular dynamics simulation showed that armchair SiCNTs can potentially
be used in separation of nitrates from water (Khataee et al., 2016).
6
SiCNTs can also be utilized in polymer composites. It has been shown that the performance of
some epoxy resins can be enhanced with SiCNTs, increasing their thermal conductivity by more
than 100% (Shen et al., 2017). Further, since SiCNTs are semiconducting, their conductivity can
be altered by molecules inside, such as HCN and CO, meaning that they could be employed in
gas detectors (Wu et al., 2008).
SiCNTs have also been studied as materials to store hydrogen (Malek and Sahimi, 2010;
Barghi et al., 2014a). Actually, Li- (Wang and Liew, 2011) and K-doped (Barghi et al., 2016a)
SiCNTs could be used to store hydrogen. However, the use of K-doped SiCNTs seemed more
eective, since such nanotubes could store up to three times more hydrogen than CNTs (Barghi
et al., 2016a). Further, SiCNTs have been used in the fabrication of membranes for reactors that
require high temperatures for catalytic membranes (Ciora et al., 2004).
Biomedical applications are also part of the potential uses of SiCNTs. SiC membranes have
been fabricated to allow the diusion of certain proteins (Rosenbloom et al., 2004). Additionally,
SiC seems to be a material compatible with living tissue (Kotzar et al., 2002).
1.4 Water conned in nanotubes
Water is a complex liquid and is considered as one of the most anomalous
uids (Pettersson
et al., 2016). It is not surprising that many mathematical models exist to try to understand
its behavior. In general, the behavior of water departs from what is expected, when compared
with other liquids (Kauzmann and Eisenberg, 1969; Debenedetti, 1996). Water density reaches
its maximum at 277 K and its viscosity decreases under pressure (Pettersson et al., 2016). The
abnormal behavior of its properties, such as isothermal compressibility, isobaric heat capacity,
and thermal expansion coecient (Kell, 1975) at temperatures between homogeneous nucleation,
231 K, and its melting point, 273 K; and the existence of at least 18 forms of ice are among the
main anomalies of water.
7
Furthermore, metastable liquid water; that is, low-density water (LDW) and high-density
water (HDW) have been identied as part of the structure in bulk supercooled water (Soper and
Ricci, 2000) and water under high pressure (Fanetti et al., 2013). The main dierences between
the two structures are in the density (about 30% dierence), and how the second coordination shell
collapses in the radial distribution function (Soper and Ricci, 2000). The transition of the two
structures can happen at temperatures below 230 K. This transition can, however, be aected by
connement in CNTs in such a way that both structures can exist at room temperature (Nomura
et al., 2017).
Under connement, hydrogen bonds can also become disturbed. Experiments have shown
that, at room temperature, fewer hydrogen bonds are formed among water molecules inside CNTs,
compared to the bulk (Ohba, 2014), and that continuous formation/breaking of the bonds occurs
at temperatures as low as 50 K (Kolesnikov et al., 2004). In polar nanotubes, such as BNNTs,
hydrogen bonds can formed between H and N, leading to a disturbance in the diusion coecient
(Won and Aluru, 2008). MD simulation show that in SiCNTs (Khademi and Sahimi, 2016) and
CNTs (Hummer et al., 2001; Cao et al., 2010) chain-like formation of water may occur, instead of
the tetrahedral structure that is common in bulk, which is due to the space limitations (Khademi
and Sahimi, 2016).
All such anomalies can be enhanced, hence changing drastically the properties of water under
connement and/or supercooling. Connement can lead to extreme phase transitions, depending
on the type and size of the nanotube: freezing can be inhibited or enhanced. Khademi and Sahimi
(2016) found that in a certain size of SiCNT water does not freeze, even at temperatures as low
as 100 K. On the other hand, Agrawal et al. (2017a) used Raman spectroscopy to nd that water
exists in an ice-like state at room temperature when conned in some CNTs.
With such disturbances to the phase transition, the questions that immediately arises are:
What is the behavior of water in CNTs and in SiCNTs, how the melting point changes, how the
8
rheological properties are aected under connement. In this work we attempt to address these
questions.
1.5 Computational Techniques and Models
In this work molecular dynamics simulation and were employed. For some specic tasks Monte
Carlo simulation was used. Without being too general, we brie
y explain important aspects of
each of these techniques in the following subsections.
1.5.1 Molecular dynamics simulation
Molecular dynamics (MD) is a computational method that allows us to study physical and
atomistic systems thorough the computation of the positions, velocities and forces that exist
among the particles or atoms in the system. The forces and interactions can be simulated by
using mathematical models called force elds. There are hundreds of models in literature for
dierent molecules and systems. Most models are developed with the purpose of simulating a
specic property of a material. For example, the TIP4P-ice (Abascal et al., 2005) water model
was developed to study ice phases and to correctly predict the freezing point of bulk water. The
same model, however, overestimates the bulk boiling point by about 50 K. Another example is the
TIP4P-2005(Abascal and Vega, 2005) water model that was developed to accurately reproduce
the maximum density of bulk water at 277 K.
Other models, such as the Terso (Terso, 1989) and the REBO (Brenner et al., 2002) poten-
tials, are \bond order" potentials, meaning that the simulated bonds between atoms are created
or destroyed dynamically, depending on the local environment. Others model are reactive force
elds, meaning that they are designed to model reactions. The RexPoN potential (Naserifar and
Goddard III, 2018) is an example of this type of potential. Some of the models mentioned above
have been used in this work and they will be explained with more detail later in this Chapter.
9
1.5.1.1 Molecular interactions and Hamiltonian mechanics
The motion of the particles in a system follows Newton's second law of motion:
F
i
m
i
B
2
q
i
Bt
2

BU
Bq
i
(1.6)
where m
i
is the mass particle, q
i
is its position, t is time and U is the potential energy, which
depends on the positions of the particles; that is, U Upq
N
q. It must be emphasized that the
potential energy is a function of the positions, q
N
q
1
; q
2
; ; q
N
only, and it is computed
based on the force eld.
To perform MD simulations, an initial system must be set up with the corresponding equations
of motions. A typical system to study is composed of N particles, and can be described by par-
ticle coordinates q
i
q
1
; q
2
; q
3
; ; q
N
and momenta p
N
p
1
; p
2
; p
3
; ; p
N
. The potential
energy Upq
N
q and the kinetic energy, Kpp
N
q
°
N
i1
|p
i
|
2
{2m
i
are used as they are part of the
Hamiltonian expression. Notice that the total energy of a system in the Hamiltonian denition is
constant. Then, the total energy of the system, in its more general form, is written as
Hpp
N
; q
N
qKpp
N
qUpq
N
qconstant; (1.7)
and the Newton's equations of motions for a particle i, in their more general form, are given by
Bp
i
Bt

BH
Bq
i
(1.8)
and
Bq
i
Bt

BH
Bp
i
(1.9)
10
We give these equations here since the statistical mechanics and the treatment of the thermostats
in MD simulations rely on such formulations. More importantly, those equations are used to
compute the positions and velocities during a MD simulation.
The systems to be simulated can be described using statistical mechanics ensembles. There are
several types of ensembles. In this work, however, there are three ensembles that are of particular
interest: The Canonical Ensemble (NVT), the isothermal-isobaric ensemble (NPT) and the Grand
Canonical ensemble (VT). The NVT ensemble is formulated such that the number of particles
(N), the volume (V) and the temperature (T) are conserved; the exact value of the total energy
is not necessarily constant. Similarly, the NPT ensemble is used then keeping the pressure P
constant. This ensemble is suitable to mimic real situations since in many experiments, the values
of T and P are known or set. These two ensembles were used in Chapters 2-4. Finally, the (VT)
uses as an input parameter the chemical potential. This is particularly useful when lling up
nanotubes with water as shown in Chapter 5.
The system with the set of 3N (due to the Cartesian nature of the positions) coupled dierential
equations can be solved numerically with a variety of algorithms. The most widely used algorithm
is the Velocity Verlet algorithm. Its popularity is due to not only its simplicity, but also its
accuracy. Moreover, the Verlet algorithm conserves the volume of the phase space of the system,
meaning that it is symplectic.
1.5.1.2 Velocity Verlet algorithm
The Velocity Verlet integration can be derived by applying a Taylor expansion on the position
vector over time. The original Verlet integration (Verlet, 1967) had a disadvantage since the
velocity and position were calculated at dierent time steps. Instead, in the Velocity Verlet
algorithm (Swope et al., 1982) they are calculated at the same time step. Thus, the Velocity
Verlet accurately calculates the velocities. The algorithm works as follows:
11
• Calculate the velocities vpt
1
2
tq based on the current positions and force from the force
eld:
vpt
1
2
tq vptq
1
2
aptqt (1.10)
• Calculate the new positions vpt tq:
vpt tq vptq vpt
1
2
tqt (1.11)
• Compute the acceleration apt tq using the force eld and the new positions qpt tq.
• Compute the velocities at the additional half time step:
vpt tq vpt
1
2
tq
1
2
apt tqt (1.12)
1.5.2 Coulombic interactions and long-range summation
The treatment of long-range interactions is a very important part of MD simulations. There are
several methods to compute the Coulombic forces. The particle-particle-particle-mesh (PPPM)
algorithm is, however, the method used in this work.
Computing the Coulombic interactions only using a typical cuto (8-12

A) may lead to artifacts
(Van Der Spoel and van Maaren, 2006) that can produce incorrect results. Nonetheless, using
a large cuto is inecient because it requires more computational time as it scales as OpN
2
q.
Hockney and Eastwood (1988) explain the use of a combination of particle-particle (using a
cuto) and a particle-mesh (using a mapping of the charges modelled as a Gaussian function in
the simulation box) that lead to a more ecient approach. The combination of the two techniques
gives rise to its name, PPPM. This methods scales as OpNlogNq.
The particle-mesh part of the algorithm treats the force as a eld quantity by positioning
the charges in a mesh. The positioning is usually done using a weight function (Hockney and
12
Eastwood, 1988). Then, the Poisson equation is solved using a discrete Fourier transform to
compute the potential eld. Moreover, the particle-particle part uses a cuto to compute the
forces among the particles within the cuto. After the forces are computed, the equations of
motions are solved to compute the new momenta and positions. A good summary of long-range
summations techniques, including the PPPM, can be found in Deserno and Holm (1998a,b).
It must be emphasized that the choice of the accuracy, how to map the charges in the mesh,
how far the charge (Gaussian base) is \spread" in the grid, and the choice of the simulation box
size are factors that help minimizing the artifacts. The key is to run several simulations with
various parameters and plot the the Coulombic energy versus the mesh size. The value that
should be taken into account should be the one when a plateau is achieved. Choosing values
where there is no constant energy can lead to artifacts.
1.5.3 Force Fields
The potentials used in this dissertation are: the TIP4P water models, the Terso potential
(Terso, 1988), the Reactive Empirical Bond Order (REBO) potential (Brenner et al., 2002) and
the RexPoN (Naserifar and Goddard III, 2018) force eld. In this subsection, a brief description
of each of them is presented.
1.5.3.1 TIP4P water models
Force elds are typically composed by two terms: bonded and non-bonded interactions. Hence,
the potential energy can be written as,
Upq
N
qU
bonded
U
nonbonded
(1.13)
13
where U
bonded
corresponds to the interactions that occurs within a molecule: bond and angle
stretching, bond and angle bending, and dihedral interactions. The U
nonbonded
accounts for the
van der Waals and Coulombic interactions among particles.
The TIP4P family only uses the non-bonded interactions, since the model simulates a rigid
molecule for computational eciency. This means that the bond and angle stretching are not
computed during a typical simulation run. These potentials are a set of 4-site rigid water models
that were developed to reproduce certain properties of water. The mathematical expression
consists of a summation of pairwise interactions of the Lennard Jones (LJ) potential (accounting
the van der Waals forces) and a Coulombic term:
UprqU
LJ
U
Coulombic
(1.14)
The LJ equations follows the typical 12-6 formulation:
U
LJ
4

r
OO

12

r
OO

6

(1.15)
Where and are the depth of the potential well and the value of r where the the potential
energy is 0. r
OO
is the distance between two O atoms.
The Coulombic term is dened as:
U
Coulombic

e
2
4
0
¸
i;j
q
i
q
j
r
ij
(1.16)
Where e is the electron charge and
0
is the permitivity of vacuum. The values of q
i
and q
j
correspond to the charge of particles i and j.
The rst TIP4P water model was developed by Jorgensen et al. (1983) using the basis of the
work done by Bernal and Fowler (1933) who was the rst to propose a 4-site water model. TIP4P
uses a rigid HOH angle of 104.52

and the length of the OH bond is xed at r(OH)= 0:9572

A.
14
The partial charges in the the H atoms are 0.52e. The negative partial charge of -1.04e is not
located on the O atom. Instead, a mass-less (M) site is located on the bisector of the HOH angle
at r(OM) = 0:15

A. This model has been widely used, specially due to its accuracy in predicting
the density 0:999 g/cm
3
at T = 25

and P = 1 atm. The experimental (expr) value is 0.9971
g/cm
3
. Later, several improvements on this model were developed while keeping the geometry
unchanged. For instance, the TIP4P-Ew (Horn et al., 2004) water model was parameterized to
be used along with Ewald summations. The boiling point (T
b
) of 370 K is accurately reproduced
(expr = 373.1 K) (Horn et al., 2005) and the diusion coecientD 2:44 x 10
9
m/s at T = 298
K and P = 1 atm. is closer to the experimental data ofD 2:07 x 10
9
m/s compared to the value
of TIP4P D 3:09 x 10
9
m/s . This was the model used in Chapter 2. Nonetheless, the above
models were not designed to model ice phases nor to match the freezing point accurately. Around
the same year that the TIP4p-Ew model was developed, the TIP4P-Ice (Abascal et al., 2005)
and the TIP4P/2005 (Abascal and Vega, 2005) water models were published. The TIP4P/2005
reproduces the diusion coecient accurately D 2:07 x 10
9
m/s at T = 298 K and P = 1
atm. However, it fails in reproducing melting and boiling temperatures, T
m
= 252.5 K (expr =
273.1) and T
b
= 401 K (expr = 373 K) respectively. Furthermore, other properties such as the
isothermal compressibility and the thermal expansion coecient are reproduced accurately. The
TIP4P-Ice model, by contrast, predicts the bulk melting point of water with high accuracy (T
m
= 272 K, expr = 273.1 K) and it produces the hexagonal ice phase in the solid phase. Although,
the viscosity reproduced of liquid water is slightly higher ( = 1.8 cP, expr = 1.0 cP at T = 298
K), the phase diagram is accurate compared to other water models. The TIP4P-Ice was used in
Chapters 3 and 4.
1.5.3.2 Bond Order Potentials: Terso and REBO
A bond order potential is an empirical interatomic model that uses the local environment
of particles to \decide" if a bond exits or not. In this models, \local environment" means the
15
neighbors of a particle and geometry that surrounds a particle. In simple words, as stated by
Terso (1988), \the more neighbors an atom has, the weaker the bond to each neighbor will be".
In bond order potentials only interactions among the bulk of the material are accounted for.
These interactions can be expressed as a sum between attractive and repulsive forces within a
cuto distance.
Since a bond order potential uses the neighbors of a particle to decide if a bond exist or not,
it necessarily has to go beyond the typical pairwise interaction. Probably this is best explained
with the following equation:
Upr
N
q
¸
i
U
1
pr
i
q
¸
i j
U
2
pr
ij
q
¸
i j k
U
3
pr
ij
;r
jk
;r
ik
q::: (1.17)
In the above equation U
1
represents interactions between one particle and external elds.
U
2
represents the interactions between 2 particles. The rst two terms are accounted in pairwise
interactions. The 3-body and higher order interactions are not part of these pairwise computation,
but they are part of many-body potentials.
One of the most famous many-body models is the Terso potential (Terso, 1988). In this
formulation the interatomic potential takes the form of (Terso, 1988):
U
¸
i
U
i

1
2
¸
ij
V
ij
(1.18)
V
ij
f
C
pr
ij
qrf
R
pr
ij
qb
ij
f
A
pr
ij
qs (1.19)
where f
R
and f
A
stand for repulsive and attractive forces (Terso, 1988). The term b
ij
is a
weight factor that depends on the environment, and the term f
C
is a cuto function that varies
from 0 to 1 as a function of the cuto r
ij
that should only include the rst shell of the neighbor
atoms (Terso, 1988). Since the environment surrounding a particle, the coordination number is
16
the most important parameter in the Terso formulation (Terso, 1988). The denitions of the
terms f
R
, f
A
and b
ij
can be found in Terso (1988).
As explained in the next Chapters, Terso potential has been widely used to model and predict
properties of solids such as Si and C. Thus, in Chapters 3 and 4 the Terso potential was used
to model SiCNTs. Nonetheless, a weakness in the potential is the incapacity of model conjugate
bonds in aromatic rings and in double bonds typical of a sp2 hybridization of pure C; similarly
it was not designed to model radicals. Brenner (1990) solved this problem by developing what is
called the Brenner potential. This potential was initially developed for hydrocarbons but has found
use in carbon nanostructures. A few years later, Brenner et al. (2002) developed the Reactive
Empirical Bond Order(REBO) potential. REBO potential is similar to the Brenner potential,
but it includes a torsion term that prevent an unrealistic rotation of the bond (Brenner et al.,
2002).
1.5.3.3 Reactive Force Field: RexPoN
Many of the bond order potentials have some capacity to model chemical reactions. This
should not be a surprise since the models have control over the existence of bonds based on the
local environment.
One of the most popular force eld to model chemical reactions is Reax, developed by
Van Duin et al. (2001) and improved by Chenoweth et al. (2008). Reax was initially con-
ceived to model chemical reactions with hydrocarbons. Inspired by the success of Reax and
by the lack of a model that can accurately predict water properties simultaneously, Naserifar
and Goddard III (2018) developed RexPoN, which is an state-of-the-art force eld based solely on
quantum mechanics with no empirical data, designed to model water. RexPoN accurately predicts
several properties of water compared to experimental data including melting point: T
m
= 273.3
K (expr = 273.15 K) and properties at 298 K: Hvap = 10.36 kcal/mol (expr = 10.52), density
= 0.9965 gr/cm
3
(expr = 0.9965), entropy = 68.4 (J/mol)/K (expr = 69.9), dielectric constant =
17
76.1 (expr = 78.4) (Naserifar and Goddard III, 2018). RexPoN includes charge polarization and
is capable of model the hydrogen bond creation/breaking process. The latest version of RexPoN
includes corrections to model a set of elements such as H, C, N, F, Cl, Br, I, P and the noble
gases (Naserifar et al., 2019). In Chapter 5, the latter version was used to model water conned
in a CNT. More details about this potential are presented in Chapter 5.
1.5.4 Monte Carlo simulation
In Chapters 2-4, the nanotubes were lled up with water by soaking them into a box full
of water. This is a completely valid procedure. Nonetheless, it might not be the most ecient
method. For this reason, using the grand-canonical ensemble combined with Monte Carlo (MC)
simulation is a good option to reduce computational time. The technique consists of an exchange
of particles between an imaginary ideal gas reservoir and the system (the nanotube), at a specic
temperatureT and chemical potential (Frenkel and Smit, 2001). The Monte Carlo part consists
of moves that can be additions, deletions, translations, or rotations of particles. The move is
accepted or rejected depending on the values of. Once the number of particles inside is constant,
the system is studied using MD simulation. This technique was used in Chapter 5.
18
Chapter 2
Complex Behavior of Ordered and Icelike Water in Carbon
Nanotubes Near its Bulk Boiling Point
2.1 Introduction
It is well-known that some
uids exhibit anomalous behavior, which is dierent from both
experimental observations and theoretical expectation for most liquids. For example, whereas
the diusion coecient of most liquids decreases under compression, there are some
uids whose
diusivity increases under the same conditions. Examples include, but not limited to, Si, (Sauer
and Borst, 1967; Kennedy and Wheeler, 1983) graphite (Togaya, 1997), and liquid metals (Cum-
mings and Stell, 1981). The most important
uid with anomalous behavior is, however, water
(Debenedetti, 2003) whose density is maximum at 277 K; it may possibly coexists within two
distinct forms (Kell, 1975) (see, however, Sellberg et al. (2014) and Laksmono et al. (2015)); it
is in a supercooled state at suciently low temperatures (Kell, 1975), and the behavior of its
isothermal compressibility, isobaric heat capacity and thermal expansion coecient (Whitby and
Quirke, 2007) at temperatures between homogeneous nucleation of 231 K and its melting point
at 273 K is unlike almost any other liquid. Moreover, it is well-known that three forms of glassy
water exist under bulk conditions, namely, low-density, high-density, and very high-density amor-
phous ice, and that, unlike most other liquids, volumetric expansion of water at 277 K and lower
19
temperatures is due to the reduction in its entropy caused by tetrahedral symmetry of the local
order around each water molecule (Kolesnikov et al., 2004) and formation of hydrogen bonds.
Given that such anomalous behavior of water occurs in an unbounded (bulk) medium, an
important question to be addressed is: what are water's properties in conned media, and in
particular in nanotubes? One cannot overstate the importance of the question, because under-
standing the properties of water in nanostructured materials and media is not only important
from a fundamental theoretical view point, but it is also critical to many physical, chemical and
biological phenomena. Experiments (Kolesnikov et al., 2004; Whitby and Quirke, 2007), as well as
computer simulations indicate that the water
ux in carbon nanotubes (Thomas and McGaughey,
2009) (CNTs) and silicon-carbide nanotubes (Khademi and Sahimi, 2011) (SiCNTs) may be three
to four orders of magnitude larger than what is predicted by the classical Hagen-Poiseuille
ow
(HP) in cylindrical tubes. In addition, it has been demonstrated that other factors, such as
functionalization of nanotubes, strongly in
uence water
ow in CNTs (Ebrahimi and Pishevar,
2015; Ramazani and Ebrahimi, 2016a). Similarly, whereas diusion of water in unbounded me-
dia follows the Einstein relation and is aected by only the local temperature and pressure, in
a conned medium the van der Waals and Coulombic interactions between water molecules and
the medium's walls in
uence their mobility, which typically reduces the rates of diusion (Zang
et al., 2009; P erez-Hern andez et al., 2012) The Stokes-Einstein relation between the diusivity
and viscosity of water also breaks down in conned media (Chen et al., 2006), including CNTs
and SiCNTs (Khademi et al., 2015; Khademi and Sahimi, 2016).
Using Raman spectroscopy, Agrawal et al. (2017a) studied the behavior of water in isolated
CNTs. They reported very high sensitivity to the diameter of CNTs and much larger increases in
the freezing transition temperature than what have been predicted theoretically. Using six isolated
CNTs of various diameters between 1.05 and 1.52 nm, Agrawal et al. (2017a) reported various
unusual phenomena that depended on the size of CNTs, including reversible melting bracketed to
105-151

C and 87-117

C, near-ambient phase transitions bracketed between 15-49

C and 3-30
20

C, and depression of the freezing point between 35

C and 10

C. Such an unusual behavior
had never been reported before.
In this Chapter we report the results of extensive molecular dynamics (MD) simulation of
water in CNTs over a wide range of temperature and with the same size that Agrawal et al.
(2017a) studied. We present the results for several important properties of water in CNTs that
provide evidence for a highly complex behavior of water, including the presence of ordered water
with connectivity much higher than its connectivity under bulk conditions, and icelike structures
that are similar to the purported ice reported by Agrawal et al. (2017a).
2.2 Computational Methods
2.2.1 Molecular Models and MD Simulation Protocol.
The MD simulations were carried out over a wide range of temperatures from 343 K to 423 K.
We built a (11,4) CNT with a C-C bond length of 1:43

A. The nanotube's length and diameter
were, respectively, ` 288

A and 10:6

A (see also below). The latter is the same as one of the
six CNTs that was used in the experiments of Agrawal et al. (2017a). The CNT was then inserted
inside a simulation cell of size 350
2
`

A
3
, lled with enough water molecules to have a water
density of 1 g/cm
3
. We then proceeded with energy minimization, followed by thermalization. The
minimization was done in three steps, using the conjugate-gradient method. First, we minimized
only the energy of the water molecules, followed by that of the CNT only, and nally for the
entire system, in order to prevent the atoms from moving large distances over short times. In
the thermalization step we increased the temperature using the Langevin thermostat, and carried
out the simulations in the (NVE) ensemble to reach the equilibrium state. The temperature
increase for all the systems that we simulated was 1 K after every 2 ps. The total time for energy
minimization and thermalization was about 1.2 ns. After reaching the desired temperature, MD
21
simulations were carried out using the velocity-Verlet algorithm in a pNVTq ensemble for 40
ns. All the simulations were carried out using the LAMMPS package (Plimpton, 1995), and the
system was visualized using Chimera (Pettersen et al., 2004). The CNT contained 3620 carbon
atoms, and there were 851 water molecules in the nanotube. The time step was always 1 fs.
The REBO potential (Brenner et al., 2002) was used to model the CNT, while the TIP4P-EW
(Horn et al., 2004) was utilized to describe the water molecules. The latter was used because
the water boiling point that it predicts, 370:3 1:9 K, is the closest to the experimental bulk
value (Horn et al., 2005). The van der Waals interactions between the carbon atoms and the
water molecules were represented by the LJ potential with a cuto distance of 11

A. The LJ
parameters for graphite (Koga et al., 2001), 3:4

A and {k
B
28:1 K
1
, were utilized
for the carbon atoms, where k
B
is the Boltzmann's factor. The Lorentz-Berthelot mixing rules,

ij
p
ii

jj
q{2 and
ij

?

ii

jj
, were utilized for the pairs ij. The Coulombic interactions
were computed using the particle-particle-particle-mesh method (Hockney and Eastwood, 1988).
2.2.2 Kirkwood g-Factor.
When simulating any phenomenon in strong connement, it is important to carry out the
computations in such a way that any possible artifact is eliminated or at least minimized. Since
we study water in nanotubes, one possible artifact is an articial crystallinity that may arise
due to the periodicity and summations methods for long-range interactions, such as the Ewald
or particle-particle-particle-mesh method (Bostick and Berkowitz, 2003). Previous studies (Van
Der Spoel and van Maaren, 2006; Yonetani, 2006; Schulz et al., 2009) showed that the Kirkwood
g-factor g
K
is an accurate tool of checking the existence of this type of artifact. g
K
is a measure
of the orientational order of the dipole moment of the molecules, and represents the correlation of
dipole
i
of water moleculei with the total dipole of the water molecules located within a distance
r from the molecule (Kirkwood, 1939). g
K
1 when the dipole moment of the molecules lacks
a specic orientation; for g
K
¡ 1 the dipole moments are parallel to an imposed eld, whereas
22
g
K
1 represents a case in which the dipole moments are aligned perpendicular to the eld
(Raicu and Feldman, 2015).
In addition to gaining insight into the problem that we study, we also used g
K
to identify the
right size of the simulation box. The g-factor is dened as (Kirkwood, 1939; Raicu and Feldman,
2015):
g
K
prq
¸
rij r

i

j

2
: (2.1)
where r
ij
is the distance between the oxygen atoms. Figure 2.1 presents the computed g
K
.
It is clear that CNTs with lengths smaller than 250

A are not appropriate. We therefore used
simulation cells of 350
2
288

A
3
.
2.2.3 The ten Wolde Parameter.
Among the many methods of analyzing solid structures, the bond-order parameters proposed
by Steinhardt et al. (1983) have been widely used. However, the problem with Steinhardt param-
eters is that they may take on nonzero values for both liquid and solid phases, making it dicult
to identify the type of structure. Due to this diculty, Rein ten Wolde et al. (1996) proposed the
use of the number of neighbors and the Steinhardt vector components to indicate the coherence
or incoherence of the \bonds" formed by the atoms of dierent molecules that, in our case, are
oxygen. Brie
y, if the dot product of the vector components from dierent molecules is 1, it
implies that they are perfectly coherently \bonded," (Rein ten Wolde et al., 1996) as is the case
in solids. However, since one may have imperfect crystals, values of 0.5 or higher are taken as
implying a solid phase (Steinhardt et al., 1983; Rein ten Wolde et al., 1996). To compute the
23
parameter we used a cuto of 3:6

A, which is the location of the rst well in the radial distribution
function. The local bond-order parameter is dened by (Rein ten Wolde et al., 1996):
q
lm
piq
1
N
b
piq
N
b
piq
¸
j1
Y
lm
r
ij
; (2.2)
whereN
b
piq is the total number of atoms, l is an integer number, m runs froml tol,Y
lm
is the
spherical harmonics, and r
ij
are the vectors between particles i and j. Since we wish to identify
possible ordered structures in the CNT, we used the global order parameterq
6
, i.e.,l 6, because
it vanishes in liquids (Steinhardt et al., 1983; Rein ten Wolde et al., 1996). Let us dene (Rein ten
Wolde et al., 1996),
d
l
pi;jq q
l
piq q

l
pjq; (2.3)
where superscript indicates complex conjugate. To every particle i one attributes a normalized
2pl 1q vector q
l
piq with components,
~ q
lm
piq
q
lm
piq

l
¸
ml
| q
lm
piq|
2

1{2
; (2.4)
Then, the histogram of values of d
l
pi;jq is used to identify the structures.
2.2.4 Mean-Square Displacement
As a measure of water diusion in CNTs, we computed the mean-square displacementsxR
2
ptqy
(MSDs) of the oxygen atoms in the axis direction:
xR
2
ptqy
1
N
N
¸
i1
rz
i
ptqz
i
p0qs
2
; (2.5)
where z
i
ptq is the axial displacement at time t.
24
2.2.5 Radial Distribution Function.
We computed the radial distribution function (RDF) gprq at the three aforementioned tem-
peratures. The RDFs were computed for the pairs O-O, H-H, and H-O pairs.
2.2.6 Mean Connectivity.
The RDF may also be used to compute the mean connectivity or coordination number Z of
the water molecules (oxygen atoms). The relation between Z and gprq is given by,
Z 4
»
rm
0
gprqr
2
dr; (2.6)
where r
m
is the rst minimum point after the largest peak in gprq, related to the smallest O-O
distance. This equation provides an estimate of the average coordination number of the oxygen
atoms, where N{V , with N being the number of the atoms, V is the volume (in

A
3
), and r
is the radial distance.
2.2.7 Cage Correlation Function.
It is well known that a property of many liquids is that, if they are cooled rapidly, they
undergo a glass transition that increases their viscosity signicantly. The viscosity increase is
not due structural changes in the liquids, but is caused by slowdown of the dynamics. The
phenomenon is usually called the cage eect, because it is as if the molecules are in a cage of
other molecules nearby. The molecules can diuse only when the cage rearranges its structure
(Gallo et al., 2010; Debenedetti and Stillinger, 2001; Bertrand et al., 2013a). To characterize the
changes in the surroundings of a molecule, cage correlation (CC) functionCptq is dened (Rabani
et al., 1997) such that, ifCptq 1, or nearly so for all times, it implies that the cage is not broken
25
and, thus, the liquid is essentially frozen. But, if Cptq decays with time, it indicates that the
cage has rearranged itself or has been broken, and that the molecules diuse signicantly in the
system, even if only anomalously slowly. The CC function has been computed by MD simulations
and measured experimentally for bulk liquids (Donati et al., 1999; Rabani et al., 1999; Kegel and
van Blaaderen, 2000). To our knowledge, Khademi et al. (2015); Khademi and Sahimi (2016)
were the rst to compute Cptq for a liquid (water) in CNTs.
We used a generalized neighbor list (Rabani et al., 1997; Weeks et al., 2000) in order to keep
track of the locations of neighbors of every atom and molecule and compute Cptq. As is well
known, a very accurate way of describing the immediate neighborhood of an atom in a water
environment is through a list of other atoms in the rst solvation shell. Then, Cptq 1 if the list
for any atom at time t remains unchanged relative to time t 0. But, if an atom or molecule
diuses outside its surrounding cage, the lists at times t and t 0 will not be identical, and
Cptq 0 at t. Thus, Cptq is computed as a function of t by averaging over all the atoms in
the CNT. We used the location of the oxygens as the basis for calculating Cptq, because it is
the closest to the center of mass of H
2
O. The aforementioned list is represented by a vector of
neighbors of each individual oxygen atoms i in the CNT with radial distance r
list
,
L
i
ptqrfpr
ij
qs; j 1; 2; ;N (2.7)
where N is the total number of atoms in the system, and
fpr
ij
q
$
'
'
&
'
'
%
1 r
ij
¤r
list
0 otherwise
(2.8)
Clearly, if r
list
is large, the computations will be long. Then, Cptq is given by
Cptq
xL
i
p0q L
i
ptqy
xL
2
i
p0qy
: (2.9)
26
To compute Cptq, we built the neighbor lists using a cuto of r
list
3:5

A, which is the lowest
value of the rst well in the RDF.
2.2.8 Intermediate Scattering Function.
The intermediate scattering function (ISF), a widely used quantity in the study of dynamics
of materials, can be measured by light or neutron scattering experiments. An advantage of
computing the ISF by MD simulation is that it can be computed for low values of the wavevector
with more accuracy, since it is usually dicult to measure the energy distribution width in
neutron scattering at such wavevectors. The ISF is dened as the spatial Fourier transform of the
van Hove function Gpr;tq,
Fp;tq
»
exppi rqGpr;tqdr: (2.10)
Since the van Hove function can be split into self G
s
p;tq (incoherent) and distinct G
d
p;tq
(coherent) parts, the same is done with the ISF:
F
s
p;tq
»
exppi rqG
s
pr;tqdr; (2.11)
F
d
p;tq
»
exppi rqG
d
pr;tqdr; (2.12)
withFp;tqF
s
F
d
, whereF
s
andF
d
are the self and distinct ISF, respectively, both of which
can be calculated directly from the MD trajectories. In particular, F
s
concerns only atoms of the
same type; it is useful for calculating the relaxation time over dierent length scales and time,
and is computed through the following equation
F
s
p;tq
1
N
C
N
¸
j1
exptirr
j
ptq r
j
p0qsu
G
(2.13)
27
where r
j
ptq is the position vector of atom j at time t, and xy implies an ensemble average.
We computed the ISF for the hydrogen atoms, as they have larger scattering area. F
s
p;tq is
typically tted to a stretched exponential function, also known as Kohlrausch-Williams-Watts
function (Kohlrausch, 1854; Williams and Watts, 1970),
F
s
p;tq exp

pt{q

; (2.14)
in which is the relaxation time.
2.3 Results and discussions
We have computed several quantities in order to characterize the state of water in the CNTs.
In what follows we present and discuss the results.
2.3.1 The Kirkwood g-Factor.
We computed the Kirkwoodg-factor for two purposes. One was to determine the appropriate
size L of the simulation cell and, hence, that of CNTs `, in order to produce results that are
independent of `. The results are shown in Figure 2.1. The second purpose was to gain insight
into the local orientational order of water that the g-factor measures. Thus, we computed the
change of the g-factor as the temperature was lowered from 423 K.
Figure 2.2 presents the results. The g-factor at high temperatures T indicates no particular
dipole orientation. AsT decreases, however, theg-factor begins to increase, indicating that there
is a dipole orientation parallel to the van der Waals force and the interaction between the CNT
and the water molecules. Because the water molecules are in a lower frequency mode (they vibrate
less) at lower temperatures, such interactions become relatively stronger, producing the increase
in the g-factor and the icelike structure. Note the gap between the g-factor for the three lowest
28
Figure 2.1: The Kirkwoodg-factor calculated using periodic boundary conditions in all directions
and the PPPM method for computing the Coulombic interactions.
Figure 2.2: The Kirkwood g-factor at several temperatures.
29
temperatures, and the higher temperatures that begin at at 373 K, indicating dierent water
structure at the lower temperatures. Note also that the highest g-factor is at 363 K. The increase
in the g-factor is consistent with the increase in the ten Wolde parameter and the slow down in
the mean-square displacements, to be described shortly.
2.3.2 The ten Wolde Parameter.
As explained in the section Methods, we computed the ten Wolde parameter d
l
pi;jq forl 6.
Figure 2.3 compares the results for water in CNTs at two temperatures with the corresponding
results for bulk water. There are very signicant dierences between the conned water order
parameter and that in the bulk. Note that values of d
6
pi;jq for about 50 percent of the water
molecules are larger than 0.50, implying that at least about 50 percent of the molecules are
coherently \bonded." In addition, at least 70 percent of the water molecules have 8-9 neighbors,
much larger than the typical 4 neighbors present in liquid water under bulk conditions. Thus, one
has a mixture of liquid and icelike states.
Figure 2.3: Comparison of the ten Wolde parameter for the bulk and conned water at T 363
K.
30
Figure 2.4: Snapshot of the distribution of water inside the carbon nanotube, with oxygen (red),
hydrogen (white), and hydrogen bonds (light blue lines).
2.3.3 Dependence of the Potential Energy on the Radial Distance.
An interesting feature is the dependence of the potential energy Uprq on the radial distance
r. As Figure 2.4 shows, there are essentially two groups of water molecules in CNT. One group
is roughly in a cylinder near the CNT's center, while the second group is near the wall. Thus,
we split the water molecules into concentric cylinders every 0:5

A in the radial direction, starting
at the center, r 0, and computed the potential energy Uprq of each group. As Figure 2.5
indicates, the water molecules near the CNT's wall possess much lower values of Uprq than those
near the center, hence indicating a slow down in the dynamics in that layer. Such a dierence in
the dynamics was also conrmed by calculating the cage correlation (CC) function, which will be
discussed next.
31
Figure 2.5: Dependence of the potential energy Uprq on the radial distance r from the CNT's
center.
2.3.4 The Cage Correlation Function.
Figure 2.6(a) presents the computed CC functions Cptq at four temperatures. Although Cptq
decays with time, but even after 40 ns it is still signicantly larger than zero. Similar to the
potential energy, we split the CNT into two cylinders. Any water molecule between r 0 (the
center) and r 1:25

A was considered as belonging to the cylinder around the tube's center. All
other water molecules belonged to what we refer to as the outer layer, which includes those near
the wall. Figure 2.6(b) presents the results for the two groups. The CC function for the outer
layer decays much slower than that near the center, hence indicating that any motion to break
the cage in the outer layer is far less likely, or much slower, than those near the center. In other
words, the dynamics only somewhat away from the center is much slower than near it. Since the
ten Wolde parameter indicates the presence of ordered structures, with 50 percent or more of the
water molecules exhibiting icelike behavior, and the molecules in the outer layer also have very
32
low potential energies than those near the center, as well as very slowly-decaying CC function, we
may conclude that an icelike phase exists in that region. Thus, two phases of water, icelike and
liquid, may co-exist in the CNT, consistent with the observation that water coexists in dierent
phases in a CNT (Barati Farimani and Aluru, 2016).
Figure 2.6: The cage correlation function Cptq. (a) The overall function at several temperatures.
(b) Comparison of Cptq for the water molecules near the tube's center and in the outer layer.
2.3.5 Intermediate Scattering Function.
Figure 2.7 presents F
s
p;tq as a function of time t and . As expected, the decay is faster
with high values of (short distances), but is slow at low wavenumbers (long distances).
Qualitatively, the data shown in Figure 2.7 are similar to those reported in several studies of
supercooled or glasslike water (Gallo et al., 1996; Chen et al., 1997; Debenedetti and Stillinger,
2001) and conned water at very low temperatures (Bertrand et al., 2013a; Diallo et al., 2015) in
various materials. The shape ofF
s
p;tq and the relaxation times that we have computed are also
comparable to those in the aforementioned studies.
As Figure 2.7 indicates,F
s
p;tq behaves qualitatively similar to the CC functionCptq in that,
it indicates a notable slowing down of the dynamics. That is, the water molecules go thorough
two stages of relaxation. The rst one is a fast process that is described by a Gaussian decay in
which the molecules move inside the cage produced by the other molecules (Gallo et al., 1996),
33
Figure 2.7: The intermediate scattering function at 343 K.
whereas the second relaxation is described by a slow stretched exponential decay, characterized
by a relaxation time
l
. This is a typical glasslike behavior where the Gaussian feature vanishes,
making the stretched exponential decay dominant as time increases, which happens for large .
In this case, at higher temperatures the function decays faster but, again, even after 40 ns it has
not vanished, except for the values that correspond to short distances. Fitting the data to Eq.
(2.14) yielded,
l
88 ns for 2, the highest peak in the structure factor.
2.3.6 Mean-Square Displacements.
The computed results for the axial mean-square displacements (MSDs) indicate a drastic slow
down in the dynamics of the water.
Indeed, if we consider the MSDs shown in Figure 2.8, we see that atT 343 K they grow with
time very slowly and after some time saturate. A similar phenomenon happens at T 353 K,
albeit after a somewhat longer time. The maximum values of the MSDs at the two temperature
34
Figure 2.8: Axial mean-square displacementsxR
2
ptqy of water (oxygen) molecules.
are around 10 12

A
2
, implying a net displacement of about 3:3 3:5

A, which is about only the
Lennard-Jones (LJ) size parameter used in various models of water, and which essentially represent
vibration of the water molecules around their positions, and not truly signicant displacements.
But, at higher temperatures, the MSDs appear to rise with time without any saturations over the
time period that we carried out our simulations.
2.3.7 Radial Distribution Function and Mean Connectivity of the Atoms.
Figure 2.9 presents the computed radial distribution functions (RDF) at T 363 K for the
pairs O-O, H-H, and H-O. The results for T 343 K and T 353 K are completely similar.
While the general trends in all the three RDFs are similar, there are also subtle dierences.
The main peak in gprq
OO
at 2:71

A is shorter than that in liquid water at room temperature,
which is 3

A. The second and the third peaks are at around 4:7

A and 5:3

A, indicating that
a central atom is surrounded by several others (whose number is computed below) in a shell of
35
Figure 2.9: The radial distribution function gprq at 363 K.
about 3:6

A, which corresponds to the rst well in gprq
OO
. gprq
OH
, on the other hand, contains
information about the hydrogen bonds and, of course, its rst peak corresponds to the O-H bond
in a water molecule, so that the second peak indicates the most likely place to nd a H atom from
an O atom. This distance is about 1:73

A. On the other hand, looking atgprq
HH
, we see that the
rst peak begins at about 1:9

A, but there are no visible wells beyond this peak. The distance is
shorter than the space between oxygen atoms, implying that the water molecules in the CNT are
very packed.
Finally, using Eq. (1) we computed the mean connectivityZ of the oxygen atoms. Z is always
between 8.7 and 9.0 in the entire range of temperatures. This is another indication of the highly
ordered structures in the CNT.
36
2.4 Conclusions
This Chapter reported the results of extensive MD simulation of water in a CNT with a
diameter 1.06 nm over the temperature range 343 K - 423 K. Various properties of water were
computed, presented, and discussed. All the computed properties indicate a strong slow down
in the dynamics of water. In the lower end of the temperature range there is strong evidence
for the presence of a two-phase system, icelike and ordered (but dynamic) water. The former
is indicated by the connectivity and the ten Wolde parameter whose values are larger than 0:5,
which is usually attributed to solid structures.
37
Chapter 3
Rheology of water in small nanotubes
3.1 Introduction
As discussed in Chapter 1, water, in addition to being vital for life, is one of the most intriguing

uids, as many of its properties exhibit anomalous behavior. For example, its density reaches
its maximum at 277 K. The anomalous behavior of several of its other properties is enhanced
as the temperature decreases (Teixeira, 2012): its isothermal compressibility increases as the
temperature decreases, reaching a maximum and then it decreases (Holten et al., 2012) at still
lower temperatures. Under supercooled condition, water density reaches a minimum at around
203 K (Mallamace et al., 2007a), while its heat capacity reaches a maximum near 250 K (Angell
et al., 2000), and its thermal conductivity attains a maximum at about 400 K (Kell, 1972). Such
intriguing anomalies at low temperatures have given rise to fundamental questions that have been
studied for at least two decades, which is why it is not surprising that numerous theoretical and
experimental studies, as well as computer simulations have been undertaken in order to understand
behavior of water over a wide range of temperature.
Compared with bulk liquid water, conned water at room temperature exhibits stronger
anomalies. For example, it is well-known that self-diusion of water strongly decreases upon
connement (Mao and Sinnott, 2000; Marti and Gordillo, 2001; Dellago et al., 2003; Kolesnikov
38
et al., 2004; Striolo, 2006; Mukherjee et al., 2007; Ma et al., 2011; Diallo, 2015; Nakaoka et al.,
2015) , and that the dynamics of the system slows down (Mamontov et al., 2006; Nanok et al.,
2009; Scal et al., 2018; Mendon ca et al., 2019). Similarly, water's thermodynamic properties may
be altered. In the past, it was reported, based on molecular dynamics (MD) simulations, that
water phase diagram can develop a \shift" when aected by connement (Kumar et al., 2005).
Later on, a phase diagram was proposed by Takaiwa et al. (2008) that indicated a decrease in the
melting point. Furthermore, dierent liquid phases, such as the low- and high-density liquid (LDL
and HDL, respectively) have been reported in carbon nanotubes (CNTs) (Nomura et al., 2017),
even at room temperatures. There is also evidence for a possible second critical point below the
freezing point, where the three phases can coexist (Chu et al., 2007; Raju et al., 2018).
Below the freezing temperature, connement may lead to inhibition of freezing and a deeply-
cooled liquid state (Bertrand et al., 2013b), giving rise to even more anomalous properties than
under the bulk conditions. Water is no exception and, indeed, analyzing the behavior of su-
percooled water in conned media below the bulk freezing temperature of 273 K reveals very
intriguing phenomena (Bertrand et al., 2013b; Diallo et al., 2015; Singh et al., 2016; Foroutan
et al., 2016). We recently showed (Khademi and Sahimi, 2011, 2016) that water in small CNTs
and their silicon-carbide counterparts (SiCNTs) of a specic size does not freeze, even below 231
K, the temperature for bulk homogeneous nucleation (Mallamace et al., 2007b). Inhibiting freez-
ing inside small tubes or pores may be accomplished, if their diameter is smaller than a critical
size of about 3-4 nm (Cobe~ na Reyes et al., 2018). Similarly, by calculating the mean-squared dis-
placements (MSDs) of the molecules in CNTs, it was determined that 8:6

A is the critical radius
in which freezing is repressed (Sakai et al., 2011). It was also demonstrated (Striolo, 2006) that
in both CNTs and SiCNTs the Stokes-Einstein relation between the viscosity and diusivity is
no longer satised, and that temperature-dependence of the self-diusivity indicates a transition
from a fragile to a strong liquid state at 230 K. Understanding the anomalous behavior of water
under such conditions is highly relevant to many chemical, biological, and physical phenomena
39
(Mashl et al., 2003), including the question of how micro-organisms survive at very low tempera-
tures, development of new instruments for preserving DNA, lubrication problems, and fabrication
of nanomaterials (Alba-Simionesco et al., 2006; Khademi et al., 2015).
Another common phenomenon in water
ow in nanotubes is the slip boundary condition on
the nanotubes' walls. In addition, Ramos-Alvarado et al. (2016) studied water-silicon interactions
and reported that the contact angle can be used to predict the slip length. Similarly, Wei et al.
(2014) reported that functionalization of conned media can weaken or strengthen the slip velocity,
which is usually orders of magnitude larger than what is predicted by the Hagen-Poiseuille (HP)
(laminar)
ow in tubes. Atomistic composition of the conned media can, in fact, drastically
change the contact angle and the slip velocity, produce very fast
ow in hydrophilic materials,
and in hydrophobic can make the layer next to the wall either to slip or stick to it; it can, of
course, also be adsorbed there (Shaat and Zheng, 2019).
The anomalous behavior of water in nanotubes has also been reported in other type of nanos-
tructured materials, such as self-assembled nanoleyers (SAMs) that are important to biological
and tribological applications. Lorenz et al. (2010), for example, found that even when the diu-
sion coecient of water increases up to two orders of magnitude when conned between layers of
alkylsilane-SAMs, no ice layers were found (although, even hydrocarbons exhibit layered struc-
tures under connement (Jabbarzadeh et al., 2006)), which is dierent from what happens in
nanotubes where ice-like behavior and layering is observed. Nonetheless, Ramin and Jabbarzadeh
(2013) reported that water does behave in an ice-like fashion when it is conned under high pres-
sures in other types of SAMs, such as n-alkanethiols, implying that despite connement of the
SAMs structure, the formation of icelike structures depends not only on the thermodynamics of
the system, but also on the atomistic composition of the nanostructured materials.
Due to the complex interactions between the water molecules, as well as between them and
a conned medium's walls, use of MD simulation is imperative. Studies involving conned water
in isolated nanotubes require, however, special treatment regarding the equilibration time. For
40
example, when the molecular structure of water is represented by such models as the SPC/E
and TIP4P/ice, the melting point cannot be reached in less than 10 ns (Fernandez et al., 2006)
of MD simulation. Moreover, the motion and arrangement of the molecules can produce false
crystalline structures, which are products of the artifacts caused by periodicity (Hub et al., 2014)
of the simulation cell, or when the Ewald summation method or the particle-particle-particle-mesh
(PPPM) approach is used to compute the electrostatic interactions, but the implementation of
the method is not done properly, which can lead to inaccuracies in the computed properties, as
demonstrated by Bostick and Berkowitz (2003).
Viscosity of
uids is usually computed by either the Green-Kubo (GK) method, or by non-
equilibrium MD (NEMD) simulation. Both are computationally expensive, due to the necessity of
requiring long simulations to generate enough
uctuating samples in the GK method, or because
of the requirement in the NEMD simulation that one must achieve steady state that is usually
reached on time a scalet that is related to the viscoelastic relaxation time of the
uid, which can
be quite long. Assuming that water is a Newtonian
uid allows use of the equations that are used
in the bulk regime or in large tubes. It was shown recently, however, that the velocity prole of
water or other Newtonian
uids
owing in nanotubes is not necessarily parabolic (the HP
ow),
even when the Reynolds number was low enough (Shaat and Zheng, 2019). This aects directly
computation of the shear rate and shear stress in nanotubes, disallowing use of the equations for
Newtonian
uids in the HP regime. Indeed, we are not aware of any study in which the rheology
of supercooled water in nanotubes, and in particular its shear stress versus shear rate diagram,
has been measured or computed. Other physical quantities, such as the intermediate scattering
function (ISF) help computing the relaxation time that is related to the viscosity (Ingebrigtsen and
Tanaka, 2018). Nonetheless, in supercooled liquids the shear-thinning behavior occurs at shear
rates several orders of magnitude smaller than the inverse of the relaxation time (Lubchenko,
2009; Furukawa, 2017). It should also be noted that the shear-thinning rheology can be observed
in any bulk liquid, as long as the applied shear rate is greater than a critical shear rate 9
c
.
41
Given the anomalous behavior of water over a wide range of temperatures, its rheology is also
important, but has remained unexplored. The goal of the present paper is to study the dynamics
and rheology of water in small nanotubes over a wide range of temperature, using extensive
MD simulations. As described below, our results indicate interesting and important rheological
behavior for water that, to our knowledge, has never been reported before.
The organization of the rest of this paper is as follows. We rst describe the molecular model
of the nanotubes and the details of the molecular simulations. The results are then presented
and discussed, including the basis for suppression of freezing of water conned in a nanotube
by computing some of the characteristic quantities. To study the rheology of conned water, we
compute the shear rate and the viscosity using a NEMD method that does not assume any specic
velocity prole for the water, and compare the results with those computed by the GK equation,
as well as with the shear viscosity of bulk water. We then discuss the dierences between the
bulk and conned critical shear rates, and their relation with the relaxation time computed by
the ISF. The paper is summarized in the last section.
3.2 Molecular Models and Computational Methods
We rst describe the nanotubes that we employed in the simulations, after which we explain
the computational procedure for the quantities have computed for the bulk and conned water.
3.2.1 The Nanotubes
We used silicon-carbide nanotube (SiCNT) as the model of conned media but, as we demon-
strated recently (Khademi and Sahimi, 2011), the same type of phenomena and results that we
present in this paper are also obtained in carbon nanotubes. It is, of course, well-known that SiC
has exceptional properties, such as thermal and mechanical stability up to very high tempera-
tures, as well as high thermal conductivity (Bai, 2011). In the past we utilized SiC for fabricating
42
nanoporous membranes (Elyassi et al., 2008) for hydrogen separation at high temperatures in
very corrosive environment, and fabricated SiCNTs as potential materials for hydrogen storage
(Malek and Sahimi, 2010; Barghi et al., 2014a, 2016a). Others have used SiC membranes for
biomedical applications (Rosenbloom et al., 2004), and studied it (Mahdizadeh and Goharshadi,
2013) as material for storage of natural gas. In the present work we used three single-wall SiCNTs,
namely, the (12,0), (20,0) and (30,0) tubes with diameters D 11:9; 19:9 and 29.8

A. Note that
among all the possible SiCNT structures, those with equal numbers of C and Si atoms were shown
(Mavrandonakis et al., 2003) to be most stable, which is also the type that we use. In setting the
nanotube's length `, one should take into account the well-known fact (Yeh and Hummer, 2004)
that the MSDs in a conned medium are in
uenced by nite-size eects. Thus, to set ` we rst
computed the MSDs of the oxygen along the axial direction z,
x
2
zptqy
1
N
N
¸
i1
rz
i
ptqz
i
p0qs
2
; (3.1)
for various lengths of the nanotube, whereN is the total number of oxygen atoms, taking into
account the change in the center of mass when computing the MSDs. While there were indeed
nite-size eects, we found that for lengths ` 300; 250 and 160

A of, respectively, the (12,0),
(20,0) and (30,0) nanotubes the nite-size eects were negligible. Thus, we set the nanotubes's
lengths at` 325; 300 and 180

A. Note that the length of the C-Si bond is 1:8

A (Menon et al.,
2004).
3.2.2 Computational protocol for bulk water
A cubic simulation cell of linear size L 35

A was lled up with 1331 water molecules,
represented by the TIP4P/ice (Abascal et al., 2005) model, which predicts a bulk freezing point
of 271.3 K. The initial density was 0.93 g/cm
3
at T 230 K and P 1 at, taken from a recent
model that predicts the density of supercooled water (Holten et al., 2012). Equilibrium MD
43
(EMD) and the NEMD simulations were carried out for, respectively, 10 ns and 50 ns. Such long
simulations are required for computing the relaxation time using the ISF by the EMD simulations,
and to attain steady state when very low shear rates are simulated in the NEMD calculations.
The NPT ensemble was used for the equilibration of the system at pressure P 1 atm. After
10 ns of equilibration the Nos e-Hoover thermostat was used in the production stage of the MD
simulations, whereas the algorithm based on the SLLOD equation of motion
1
was employed in
the NEMD (Evans and Morriss, 1984). The self-part of the ISF was computed using
F
s
p;tq
1
N
C
N
¸
j1
exptirr
j
ptq r
j
p0qsu
G
; (3.2)
where r
j
ptq is the position vector of atom j at time t, and xy implies an ensemble average.
We computed the ISF for the hydrogen atoms, as they have larger scattering area. F
s
p;tq is
typically tted to a stretched exponential function, also known as Kohlrausch-Williams-Watts
function (Kohlrausch, 1854; Williams and Watts, 1970),
F
s
p;tqa exp

pt{q

; (3.3)
in which is the relaxation time, and a and are constants. In the NEMD simulation the
shear viscosity was computed using

xy
9
; (3.4)
where
xy
is the shear stress, and 9
is the shear rate. We used a wide range of the shear rates
in order to identify the onset of non-Newtonian behavior of water, if any.
1
The name SLLOD is due to the use of the transposed Doll's tensor, the dyadic product of the positions and
momenta, which was named after the Kewpie Doll by Wm. G. Hoover.
44
3.2.3 Computational protocol for conned water
After inserting the nanotubes in individual simulation boxes with linear dimensions of 400, 350
and 240

A for the (12,0), (20,0), and (30,0) tubes, respectively, water molecules were distributed
in the box, including in the nanotube, with a density of 1 g/cm.
3
The SiCNTs were modeled
by the Terso (Terso, 1989) potential that has been widely used for estimating the mechanical
properties of SiCNTs (Zhang and Huang, 2008; Setoodeh et al., 2009; Shen, 2009; Memarian et al.,
2017). The interactions between water and the SiCNTs were represented by the Lennard-Jones
(LJ) potential, whose parameters for the SiCNTs were taken from Malek and Sahimi (Malek
and Sahimi, 2010), while the partial charges were taken from the ab initio calculations reported
previously (Mavrandonakis et al., 2003). The Lorentz-Berthelot mixing rules were employed to
compute the LJ parameters of pairs of atoms.
The energy of the system was then minimized, which typically took between 1.5 and 3 ns.
Then, EMD simulations were carried out at each temperature in the NPT ensemble at P 1
atm for 10 ns in order to allow the water molecules to move in and out of the nanotubes. After
equilibrium was reached and the potential energy was stable, water molecules outside the nanotube
were removed, the axis of the simulation box was aligned with that of the nanotubes, and periodic
boundary conditions were imposed. The system was then taken to its nal temperature. We also
ran a set of simulations in the (12,0) SiCNT at 80 K and 290 K in order to visualize solid and
liquid states on both sides of the freezing point.
After the structures were generated, separate sets of MD and NEMD simulations were carried
out. The EMD simulations were done for 10 ns to compute the CCF, the MSDs, and the radial
distribution function (RDF) gprq for at several temperatures. Calculation of the ISF, however,
required 200 ns to reach its low values. In order to classify the molecules as high-, low-, and
normal density water (HDW, LDW, and NDW, respectively), we used the Voronoi tessellation
method to compute the volume occupied by each type. The local density of the water molecules
45
is given by the inverse of the volume of the Voronoi blocks. Thus, we classied water as follows:
water molecules with a density larger than 1.0 g/cm
3
were considered as the HDW; molecules
with density between 0.9 and 1.0 g/cm
3
were labeled as the NDW, whereas the molecules with
a density below 0.9 g/cm
3
were treated as the LDW. We then computed gprq for the LDW and
HDW.
To compute Cptq, the CCF, we used a generalized neighbor list (Rabani et al., 1999) in order
to keep track of the locations of neighbors of every molecule, using the rst solvation shell as the
boundary, since we wish to understand how fast the cage rearranges itself. If Cptq 1, the list
for any atom at time t remains unchanged relative to time t 0. If, however, a molecule diuses
outside its surrounding cage, the lists at times t and t 0 will not be identical, and Cptq 0
at time t. Thus, Cptq was computed by averaging over all the atoms in the nanotube. We used
the location of the oxygen atoms as the basis for calculating Cptq, because it is the closest to the
center of mass of H
2
O. The aforementioned list was represented by a vector of neighbors of each
individual oxygen atoms i in the CNT with radial distance r
list
,
L
i
ptqrfpr
ij
qs; j 1; 2; ;N ; (3.5)
where N is the total number of atoms in the system, and
fpr
ij
q
$
'
'
&
'
'
%
1 r
ij
¤r
list
0 otherwise:
(3.6)
We used a cutor
list
3:25

A, which is the lowest value of the rst well in the RDF. Using a
larger value of r
list
may lead to misinterpreting the results. Then, Cptq is given by
Cptq
xL
i
p0q L
i
ptqy
xL
2
i
p0qy
: (3.7)
46
The MSDs in the axial direction z were computed through Eq. (1). The zero-shear rate
viscosity was computed using the GK equation, by integrating the stress-stress autocorrelation
function of the three o-diagonal entries of the viscous pressure tensor:
pTq
V
3k
B
T
»
8
0

x
xy
ptq
xy
p0qy
x
xz
ptq
xz
p0qyx
yz
ptq
yz
p0qy

dt;
(3.8)
where V is the volume of the system, and k
B
is the Boltzmann's constant.
The NEMD simulations were carried out for between 10 and 20 ns. At time t 0 an axial
force was applied to every water molecule in order to induce motion. To do so, the zaxis (axial)
degree of freedom was removed from the thermostat to avoid articial heating of the system. It
took about another 8 ns for the systems to reach steady state. Once the velocity of the center of
mass of water molecules reached steady state, the external applied force was turned o to allow
water molecules to decelerate. Computation of the deceleration could have been done by various
ways. Here, however, we only used up the 800 fs right after the external force was removed to avoid
inaccuracies due the motion of the zaxis degree of freedom of the thermostat. Such approach
has been used before in several studies of conned water (Chen et al., 2008; Xu et al., 2011; Liu
et al., 2012; Kheirabadi et al., 2014; Liu et al., 2016a,b) Data for the velocities were collected
every 1 fs, and averaged over 5 fs intervals, where the average was dened by, vM
1
°
N
i
m
i
v
i
,
withM being the total water mass in the tube,m
i
andv
i
the mass and velocity of water molecule
i, and N the total number of water molecules.
47
To compute the quantities of interest, we proceeded as follows. After removing the external
force, the deceleration a of the water molecules was computed. Then, the shear stress at the wall
was calculated by

M a
DL
; (3.9)
which is obtained by a force balance and is independent of the rheology. Note that because the
LJ cuto was set at 12, 20 and 30

A for the three nanotubes, which are slightly larger than each
nanotube's diameter, all water molecules interact with the walls. Thus, Eq. (9) is representative
of the shear stress in the nanotube. We also point out that Eq. (9) does not contain a pressure
term, because the force applied to each particle to induce
ow was used. To compute the shear
rate 9
, we divided the interior of the nanotube into ve concentric cylinders, the maximum number
of layers in which the water molecules were distributed in the tube. The center-of-mass velocity
for the groups was used for computing the mean shear rate, given by

C
BV
z
Br

BV
r
Bz
G
: (3.10)
Thus, sinceBV
r
{Bz 0, we have

B
BV
z
Br
F

B
V
z
r
F
; (3.11)
where V
z
{r denotes the dierence in the values between two layers at a distance z.
Finally, the average shear rate was computed using those obtained by Eq. (11). Multiple values
of the applied force were employed to simulate water
ow in the nanotubes over a wide range of
the
uid velocities.
Let us point out that, due to slip on the nanotube's wall, computing the shear rate is a dicult
problem, and particularly so if its value is small. For example, it was not possible to compute
48
the shear rate at velocities v
z
smaller than 10 m/s, which are velocities that are comparable
with the thermal velocities of the molecules. Moreover, the water layers in the nanotube may
have velocities that are close to each other. To overcome this diculty, the velocity of the center
of mass of each layer was averaged over all the time between the beginning of the steady-state
velocity prole and the time at which the external force was removed. This removed most of
the thermal noise, making it possible to obtain accurate results. In addition, a minimum of 30
realizations were generated for each data point to improve the statistics, and for velocities smaller
than v
z
30 m/s at least 50 realizations were used.
In both the EMD and NEMD simulations, long-range Coulombic interactions were computed
by the PPPM method (Hockney and Eastwood, 1988). The SHAKE algorithm (Ryckaert et al.,
1977) was used to keep the bonds and angles in the water molecules intact. Temperature was
adjusted by the Langevin thermostat, and its increase to higher values, when needed, was done at
a rate of 1 K after every 5 ps. All the simulations were carried out using the LAMMPS package
(Plimpton, 1995), and the systems were visualized using the VMD package (Humphrey et al.,
1996). The time step was always 1 fs. As mentioned earlier, water molecules were represented by
the TIP4P/ice (Abascal et al., 2005) model. Their viscosity at 273 K, calculated by the TIP4P/ice
model, is 1.78 cP Louden and Gezelter (2018), which is larger than the experimental value of 1
cP, but is of the correct order of magnitude.
3.3 Results and discussion
We computed several important physical quantities that characterize the rheology and dynam-
ics of supercooled water in small nanotubes. In what follows we present and discuss them.
49
3.3.1 The cage correlation function and mean-squared displacements
Figure 3.1 presents the CCF at three widely disparate temperatures. At 80 K the CCF does
not change, indicating a frozen state. It decays clearly at 230 K, indicating that the cage formed
by the water molecules rearranges itself continuously due to jumping of the molecules into a
dierent cage. The decay of the CCF is the sharpest at 290 K, as expected.
10
6
10
5
10
4
10
3
10
2
10
1
10
0
10
1
time (ns)
0.0
0.2
0.4
0.6
0.8
1.0
C(t)
T (K)
80
230
250
270
273
290
Figure 3.1: The cage correlation function Cptq. The results at 80 K indicate a frozen state,
whereas those at other temperatures are representative of a liquid state.
Figure 3.2 presents the results for the MSDs. They are practically zero at 80 K, not chang-
ing over time. But, they grow continuously with the time at 230 K and higher temperatures.
Moreover, at longer times the MSDs vary essentially linearly with the time, indicating ordinary
diusive motion. The inset in Fig. 2 shows the diusion coecient D, computed based on the
50
MSDs. The bulk values were taken from Weiss et al. (2011). The results in Figs. 1 and 2 are
completely consistent with the existence of a liquid state below the bulk freezing temperature.
0 2 4 6 8 10
time (ns)
0
25
50
75
100
125
150
175
200
<
2
z(t)> (Å
2
)
T (K)
80
230
250
270
273
290
220 240 260 280 300
T (K)
10
10
10
8
10
6
10
4
D (cm
2
/s)
T (K)
Confined
Bulk
Figure 3.2: Axial mean-squared displacementsx
2
zptqy of water. The results at 80 K indicate a
frozen state. The inset compares the diusivities in connement with the bulk values. The bulk
values were taken from Weiss et al. (2011)
3.3.2 Radial distribution function and density
uctuations of liquid
water
The radial distribution function gprq of the oxygen-oxygen pairs is shown in Fig. 3.3(a). The
dierences in the relative heights of each peak atT 80 and 230 K are consistent with the results
shown in Figs. 3.1 and 3.2. The second and third peak at 230 K are more \
at," when compared
with those at 80 K, indicating that at 230 K one has a liquid state.
51
The presence of the three types of density was conrmed by the computation using the Voronoi
volumes, as described earlier, as well the RDF gprq. Locally low and high densities have distinct
gprq patterns. This is shown in Figs. 3.3. For high density the second peak collapses, with the
rst one leaving a small shoulder at around r 3

A. Moreover, the rst peak corresponding to
the high density develops before the rst peak representing the low density, which is due to the
closeness of the water molecules in the former region. The second peak corresponding to the low
densities is broad and larger than the corresponding one in the region with high density, due to
the larger probability of nding water molecules than in the latter.
At 230 K the normal density water tends to occupy the largest portion among the three,
representing about 55% of the total volume, followed by the high and low densities at, respectively,
33 and 12 %. The fractions change only slightly at 290 K, with the values for the regions with
normal, high, and low densities being, respectively, 59, 31, and 10 %, with all the changes being
well in the range of the numerical error. A summary of the density
uctuations in the nanotube
is shown in Figs. 3.4. At T 230 K and 290 K, the low density represents the lowest portion,
and there is a peak at around a density of about 0.925 g/cm,
3
which closely matches the density
of bulk water at the same temperature (Holten et al., 2012). In contrast, at T 290 K there
is not a unique peak. Instead, the molecules with normal density have a broader distribution.
Nevertheless, the fractions of the three types of local densities at 290 K in the larger nanotubes are
signicantly dierent. Higher volume fractions of the normal density, 63 % in the (20,0) and 70
% in the (30,0) nanotubes, are present. In contrast, the fractions of the regions with high density
in the (20,0) and (30,0) nanotubes are only 28 and 23 %. The low density type also experiences
such a decrease, which are 9 and 7 %. The anisotropy produced by the connement can be a
possible explanation of the rheology, as explained below. Note that structural changes in sheared
liquids have been interpreted as the ngerprint of non-Newtonian rheology (Jabbarzadeh et al.,
2006; Ingebrigtsen and Tanaka, 2018).
52
0.0 1.0 2.0 3.0 4.0 5.0 6.0
0
1
2
3
4
5
6
7
8
T (K)
80
230
273
1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5
0
1
2
3
4
5
g(r)
T (K)
230
250
270
290
1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5
r (Å)
0
1
2
3
4
5
T (K)
230
250
270
290
(a)
(b)
(c)
Figure 3.3: The radial distribution function gprq, which provide evidence for a (a) liquid state,
and the existence of (b) low- and (c) high-density water regions in the (12,0) nanotube
53
0.80 0.85 0.90 0.95 1.00 1.05 1.10
290 K
0.80 0.85 0.90 0.95 1.00 1.05 1.10
Density (g/cm
3
)
0.00
0.02
0.04
0.06
0.08
0.10
0.12
0.14
0.16
0.18
Relative frequency
230 K
Figure 3.4: Histogram of density
uctuations of water in the (12,0) nanotube
3.3.3 Intermediate scattering function and the relaxation time
Figure 3.5 presents the ISF for bulk and conned waters. Note that in the MD simulations
of water at low temperatures that we consider, the bulk water does not freeze, because in order
to produce ice one must add nucleation sites for ice; otherwise, it is very unlikely that MD
simulation of bulk water at such temperatures can produce ice. In other words, one has to put
the supercooled water in contact with, for example, hexagonal ice in order to produce the frozen
state; see for instance, Naserifar and Goddard III (2019), who recently simulated supercooled bulk
water without ice.
The relaxation times that were computed for bulk water at 230 K, 250 K, 270 K and 290 K
were, respectively, 60 ps, 23 ps, 12.5, and 6.9 ps. The corresponding values for conned water
in the (12,0) nanotube at the same temperatures were 4.18 ns, 2.08 ns, 1.01 ns and 0.57 ns,
respectively, about two orders of magnitude larger than those of bulk water. The fact that the
relaxation times of conned water are much larger than those in the bulk provides evidence for
the existence of a possible critical shear rate 9
c
for non-Newtonian rheology. As mentioned earlier,
one usually nds that 9
c
for supercooled liquids is smaller than the inverse of the relaxation time,
54
implying that 9
c
is very small. Here, we present evidence indicating that the critical shear rate of
conned water is smaller than the corresponding value for the bulk water, implying that shear-
thinning rheology in connement sets in much earlier than in the bulk. We shall return to this
point shortly.
10
5
10
3
10
1
10
1
10
3
T (K)
230
250
270
290
(20,0) 290
(30,0) 290
10
6
10
5
10
4
10
3
10
2
10
1
10
0
10
1
time (ns)
0.0
0.2
0.4
0.6
0.8
1.0
F
s
( , t)
T (K)
230
250
270
290
(a) (b)
Figure 3.5: The intermediate scattering functionF
s
p;tq for (a) bulk water, and (b) water conned
in the nanotube.
3.3.4 Rheology
The viscosity of bulk supercooled water was computed by the NEMD simulation. The results
are shown in Fig. 3.6. We note that the onset of non-Newtonian rheology occurs for shear
rates between 1 and 10 ns,
1
which is about one order of magnitude lower than the inverse of
the relaxation times. The eective zero shear-rate viscosity of the water conned in the (12,0)
nanotube was computed by the GK formalism. The results are shown in Fig. 3.7, and indicate
that water in the nanotube is a very viscous
uid at 230 K. The viscosity decreases gradually as
temperature rises.
Zero shear-rate viscosity of water at room temperature, computed by the GK relation, was
also reported previously by us (Khademi and Sahimi, 2011). Our estimates are about twice larger
55
10
8
10
9
10
10
10
11
10
12
(s
1
)
0
5
10
15
20
(cP)
T (K)
230
250
270
290
Figure 3.6: Viscosity of bulk water as a function of the shear rate, computed by the NEMD
simulation.
230 240 250 260 270 280 290
T (K)
5
10
15
20
25
30
35
40
45
(cP)
Figure 3.7: Zero-shear rate viscosity of water in the nanotube, computed by the Green-Kubo
equation.
than the experimental values reported for bulk supercooled water (Dehaoui et al., 2015) over the
same temperature range. Note that, as has also been pointed out by others (Xu et al., 2011;
56
Liu et al., 2012), use of the GK relation in conned media is not without problems, because the
equation was originally derived for homogeneous
uids. Moreover, the frequency of sampling and
the periodicity may generate artifacts. Despite this, as explained below, the calculations based
on the GK relation, can still be useful.
To obtain more accurate estimates of the viscosity based on sound theoretical foundations,
we also carried out a separate set of NEMD simulations to compute shear-rate dependence of
the eective viscosity. Figure 3.8 presents the shear stress-shear rate diagram for water in the
nanotubes as a function of the temperature. In all cases the shear stress does not vary signicantly
with the shear rate 9
when 9
is small. But, as 9
increases, so also does the shear stress. The
results shown in Fig. 3.8 indicate clearly a non-Newtonian rheology at all temperatures.
10
5
10
6
10
7
10
8
10
9
10
10
(s
1
)
10
4
10
3
10
2
10
1
10
0
10
1
(MPa)
T (K)
230
250
270
290
(20,0) 290
(30,0) 290
Figure 3.8: Temperature- and shear rate-dependence of the shear stress in the nanotube, computed
by the NEMD simulation.
To understand this better, as well as obtaining more accurate estimates of water viscosity in the
nanotube, we dene an eective viscosity by Eq. (4). For Newtonian
uids, is independent
of the shear rate and depends only on T (and pressure P ) and the molecular structure of the
solution, but for non-Newtonian
uids depends also on 9
. The results for, presented in Fig. 3.9,
57
provide further evidence that the connement and the supercooled state of water enhance the non-
Newtonian rheology, and that water resembles a shear-thinning
uid since its viscosity decreases
with the shear rate 9
. Note also that as the temperature and the nanotubes' diameter increase,
the dependence of on 9
weakens, hence indicating that rheology of water in the nanotube at
temperatures higher than room temperature resembles increasingly that of Newtonian
uids.
We note that the even though it is not possible to precisely determine the critical shear rate
for conned water, it is evident from Fig. 3.9 that it is at least two or three orders of magnitude
smaller than the inverse of the relaxation times in the same systems. This indicates that shear-
thinning rheology in nanotubes emerges earlier than in the bulk.
10
5
10
6
10
7
10
8
10
9
10
10
(s
1
)
0.0
2.5
5.0
7.5
10.0
12.5
15.0
17.5
(cP)
T (K)
230
250
270
290
(20,0) 290
(30,0) 290
Figure 3.9: Temperature- and shear-rate dependence of the eective viscosity of water, computed
by the NEMD simulation. They indicate non-Newtonian rheology.
The eective viscosity for the smallest (zero) shear rate that we computed by the NEMD
simulations in the nanotube are smaller than that computed by the GK formalism, but much
closer to the experimental data (Dehaoui et al., 2015) for supercooled water under bulk condition
over the same temperature range. We explain this by noting that (Zaragoza et al., 2019) thermal

uctuations of shear stress at equilibrium will not only be aected by the viscosity, but also by
58
the slip boundary condition at the wall, which is known to exist in
ow through small nanotubes.
Thus, because the GK equation was intended for homogeneous
uids, it does not take into account
the slip on the walls.
Therefore, the GK provides overestimates of the eective viscosity in conned media. This
implies that either the GK formulation must be modied in order to take into account the slip
eect, or that the estimates that it provides for the eective viscosity of liquids in conned media
with slip boundary conditions should be viewed as a sort of upper bound to the true values, as it
ignores the slip on the walls. Indeed, recent research (de la Torre et al., 2019) has shown that for
conned
uids between parallel plates the GK equation is inadequate for computing the viscosity.
Finally, we note that although the viscosities for the bulk and conned water were computed
by dierent methods, the dierences between the relaxation times agree with the slow dynamics
encountered due to connement, which was reported previously in several studies (Kell, 1972;
Mao and Sinnott, 2000; Marti and Gordillo, 2001; Dellago et al., 2003; Kolesnikov et al., 2004;
Striolo, 2006; Kumar et al., 2005).
3.4 Conclusions
This chapter reported on a study of rheology of supercooled water in small nanotubes. We
reported the results of extensive equilibrium and non-equilibrium molecular dynamics simula-
tions of supercooled water in small SiCNTs over a wide range of temperatures. The properties
computed, namely, the CCF, MSDs, and gprq, shear stress-shear rate diagram, and the eective
viscosity, indicate that not only does water not freeze at temperatures way below the bulk freezing
temperature of 273 K, but also that regions with high and low densities coexist in the nanotube at
such temperatures and above the freezing point. The volume fractions of the regions with locally
high and low densities tend, however, to decrease as temperature and diameter of the nanotube
increase. We also note that the non-Newtonian regime tends to be in a narrower range of shear
59
rates as the nanotube's diameter increases. Moreover, while the presence of the liquid with dif-
ferent densities in smaller nanotubes is very clear, the same is not true in larger ones, where the
non-Newtonian behavior is weaker.
We also reported clear evidence for non-Newtonian rheology of conned water over the same
range of temperature at high velocities. We suspect that the existence and the relative abundance
of the water with dierent local densities may be linked to the non-Newtonian rheology, since
it has been shown that the structure of liquids has a strong relationship with their viscosity,
especially in supercooled water (Ingebrigtsen and Tanaka, 2018). Most recently, it was proposed
by Naserifar and Goddard (2019) that water could behave as a polymer chain over certain time
scales, which may provide further evidence for its non-Newtonian rheology in connement.
These results may have important implications for
ow of water and similar liquids in nanos-
tructured materials, for biological systems, and for the design of devices for various medical
applications at low temperatures.
60
Chapter 4
Universal Intrinsic Dynamics and Freezing of Water in
Small Nanotubes
4.1 Introduction
Static and dynamic properties of conned water dier signicantly from those in the bulk
(Mondal and Bagchi, 2019; Cobe~ na Reyes et al., 2018; Scal et al., 2018; Tomo et al., 2018; K ohler
and Bordin, 2018; Chakraborty et al., 2017; Ramazani and Ebrahimi, 2016a; Chong and Ham,
2016; Barati Farimani and Aluru, 2016; Khademi and Sahimi, 2016; Ramazani and Ebrahimi,
2016b; He et al., 2013; Hernandez-Rojas et al., 2012), which have opened up the possibility of
exploiting them for a large number of potential applications. In particular, properties of water
in nanotubes have been studied extensively because, for example, water-nanotube systems may
be used for desalination processes (Dong et al., 2018; Ebrahimi et al., 2018; Tunuguntla et al.,
2017; Gao et al., 2017; Bhadra et al., 2016; Tu et al., 2016; Zhao and Wu, 2015; Das et al., 2014;
Azamat et al., 2014; Hughes et al., 2012; Mishra and Ramaprabhu, 2010), in wastewater treatment
(Dokmaj et al., 2020; Eskhan et al., 2019; Ou et al., 2019; Mubarak et al., 2014), drug delivery
(Kaur et al., 2019; Kafshgari et al., 2019; Larnaudie et al., 2018; Karimi et al., 2015; Martincic and
Tobias, 2015), and many other important applications (Ebrahimi et al., 2020; Shimizu et al., 2019;
61
Mahbubul et al., 2018; Ge et al., 2017; Abdullah and Kamarudin, 2017; Sahimi and Ebrahimi,
2019).
Connement slows down the dynamics of water conned in nanostructured materials, such
as carbon nanotubes (CNTs) (Cobe~ na Reyes et al., 2018; K ohler and Bordin, 2018; Chakraborty
et al., 2017; Chong and Ham, 2016; Varanasi et al., 2018; Pugliese et al., 2017; Agrawal et al.,
2017a; Diallo et al., 2015), as well as silicon-carbide nanotubes (Cobe~ na-Reyes and Sahimi, 2020a;
Khademi et al., 2015; Farmahini et al., 2014; Taghavi et al., 2013; Yang et al., 2011; Khademi and
Sahimi, 2011) (SiCNTs) and boron-nitride counterpart (Won and Aluru, 2008). For example, the
relaxation time of conned water increases by several orders of magnitude (Ebrahimi et al., 2020;
Sahimi and Ebrahimi, 2019; Cao et al., 2018; Briganti et al., 2017; Pajzderska et al., 2014), while its
diusivity decreases (Cobe~ na Reyes et al., 2018; K ohler and Bordin, 2018; Zheng et al., 2012). The
increase in the relaxation times (Cobe~ na-Reyes and Sahimi, 2020a), disruption in the hydrogen
bond (HB) network (de Almeida Ribeiro and de Koning, 2020), and short lifetime of polymeric
behavior of water (Naserifar and Goddard, 2019) have also been proposed as the reasons for water
to behave as a non-Newtonian
uid when supercooled (de Almeida Ribeiro and de Koning, 2020)
or under connement (Cobe~ na-Reyes and Sahimi, 2020a). Another outstanding characteristic of
water in small nanotubes is the break down of the Stokes-Einstein relation between diusivity
and viscosity in supercooled (Ansari et al., 2020; Dubey et al., 2019; Shi et al., 2018; Kawasaki
and Kim, 2017; Xu et al., 2009) and conned water (Zaragoza et al., 2019; K ohler et al., 2017;
Khademi and Sahimi, 2011; Mallamace et al., 2010), an anomaly in the diusion coecient that
implies the existence of two distinct diusion regimes, separated by a crossover temperature that
is approximately 230-240 K at which the fragile-to-strong (FTS) transition occurs (Dubey et al.,
2019; Shi et al., 2018; Khademi and Sahimi, 2011; Mallamace et al., 2010).
Water can, of course, freeze and form structures with distinct shapes. In particular, ice-like
structures form in CNTs with a number of dierent shapes (Takaiwa et al., 2008). Using molecular
dynamics (MD) simulations, Takaiwa et al. (2008) reported formation of nine distinct ice types in
62
several CNTs. The highest melting point of such ices was about 290 K, formed by water conned
in a CNT of diameterD 10:8

A, whereas the lowest melting point was around 180 K in a CNT
with a diameter of 9:8

A (Takaiwa et al., 2008). Similarly, Raju et al. (2018) studied formation of
one- and two-dimensional ice in parallel graphene sheets and CNTs using the ReaxFF force eld
(Raju et al., 2018). They reported a rst-order transition whereby hexagonal ice was present in
one of the CNTs, as well in one of the water-graphene systems (Raju et al., 2018). But they also
claimed that a continuous phase transition occurs based on the behavior of the radial distribution
function and the changes in the diusion coecient, and that the continuous phase transition
indicated the existence of a supercritical system (Raju et al., 2018). On the other hand, K ohler
and Gavazzoni (2019) reported distinct behaviors of water conned in MoS
2
nanotubes. In the
smallest nanotube with diameter of D 9

A water had a higher mobility than in a nanotube of
diameter D 12:5

A, where frozen water was present.
Nonetheless, experiments have also provided evidence for changes in the melting point of
conned water. Using Raman spectroscopy, Agrawal et al. (2017a) reported extreme sensitivity
of the freezing point of water conned in CNTs to their diameters. The melting point in a CNT
of diameter D 10:5

A was reported to be between 378 K and 424 K, whereas for a CNT with
D 11:5

A the possible range was between 238 K and 283 K, and when D 12:4

A, the melting
point was below 193 K. These results were supported by MD simulations (Cobe~ na Reyes et al.,
2018; Pugliese et al., 2017). A more recent paper by Chiashi et al. (2019) indicated that the
changes in the melting point might not be as drastic as reported by Agrawal et al. (2017b). Using
photoluminescence spectroscopy, Chiashi et al. (2019) reported that in most of the CNTs that
they experimented with the melting points were close to that in the bulk. They also reported that,
depending on the diameter of the nanotube, the same type of ice may have a dierent melting
point.
The search for the change in the melting point of water in nanostructured materials, by either
experiment or MD simulation, is a dicult problem. Depending on the material and diameter of
63
the nanotube, complete suppression of freezing might occur or, in some cases, it can be dicult
to determine the freezing point, especially if it is too low. Recent MD simulations (Khademi and
Sahimi, 2016) indicated that water may still behave as a liquid in SiCNTs at temperatures as low
as 100 K. Similarly, Agrawal et al. (2017a) did not observe a freezing transition in an isolated
CNT of diameter D 12:4

A at temperatures in the range 190 K - 390 K, and Chiashi et al.
(2019) did not report melting in one of their CNTs. This implies that at least in some cases the
melting point may be shifted well below the bulk freezing point.
Connement also aects the intrinsic dynamics of the HB network and the orientation of the
dipoles of the water molecules. For example, the number of HBs can be altered by the type of
the atoms forming the nanotube, or by ionic particles in contact with water (Wang et al., 2020).
The HBs can form and break faster, depending on the environment, or they can be broken upon
deformation of the nanotubes' walls (Mendon ca et al., 2019). Moreover, the number of the HBs
may decrease, depending on the charge and interaction parameters of the atoms that form the
nanostructures (Zhang et al., 2017; Hernandez-Rojas et al., 2012). In addition, there seems to
be a relationship between the dipole orientation and the HB network (Ma et al., 2017), where
the typical dipole orientation is measured with respect the axial direction, z, which in CNTs is
about (Foroutan et al., 2020; Zeng et al., 2018) 30

. Nonetheless, can be altered, depending
on the nanotube and the partial charges of its constituent atoms. Mei et al. (Mei et al., 2015)
reported that not only is the dipole angle of water conned in graphite layers dierent from 30

when the partial charges of the carbon atoms are [articially] altered, but that water can form
fewer HBs compared to the case in which the charges in the graphite walls are not altered. After
all, conned water will form structures that maximize the number of HBs, but if they are weak,
the melting and boiling points would be depressed (Takaiwa et al., 2008).
Given the above discussion, the questions that arises are whether the changes in the melting
temperature are universal, independent of the type of nanotube, and how the dynamics of water
is related to such changes. In this paper we address these questions by using extensive MD
64
simulations of water conned in SiCNTs and CNTs with diameters similar to those of the CNTs
used by Chiashi et al. To study the dynamics and determine the melting point, we have computed
several characteristic quantities.
This chapter is organized as follows: First, we describe the computational details of the MD
simulations. Next, the results are presented and discussed. The chapter is summarized in the last
section.
4.2 Molecular models and details of the computations
4.2.1 Computational models
Two types of nanotubes, namely, CNT and SiCNT, were employed in our MD simulations,
although most of our calculations were with the latter ones. The extraordinary properties of SiC,
such as their thermal and mechanical stability (Pokropivnyi and Silenko, 2006; Bai, 2011) up to
very high temperatures, as well as their resistance to corrosive environments (Barghi et al., 2014b)
are, of course, well known. We have also fabricated SiCNTs for hydrogen storage (Barghi et al.,
2014b, 2016b), gas adsorption (Malek and Sahimi, 2010) and separation (Elyassi et al., 2007).
Others have studied SiC as material for storage of natural gas (Mahdizadeh and Goharshadi,
2013) and for electronic applications (Zhang et al., 2019; Talla, 2019). More recently, SiCNTs
have exhibited compatibility with phospholipd membranes (Raczyn ski et al., 2019), implying that
they may be used for carrying drugs to biological cells (Mehrjouei et al., 2017).
Six single-wall SiCNTs of various diameters, lengths and chiralities, namely, (9,2), (6,6), (8,4),
(10,2), (10,3) and (11,2), were used with diameters that were, respectively,D 10:07, 10.31, 10.5,
11.05, 11.70 and 12:3

A, matching closely those of the CNTs used by Chiashi et al. The lengths`
of the SiCNTs were, in the same order as above, 330, 328, 322, 361, 382, and 350

A. Two CNTs,
(10,5) and (11,3), with diameters D 10 and 10:3

A and lengths ` 326 and 338

A were also
65
used in the MD simulations in order to compare the decay of the HB correlation function (HBCF)
in the two types of nanotubes, as described below.
To ll up the nanotubes with water we placed the (9,2), (6,6), (8,4) SiCNTs and the two CNTs
in individual boxes of linear dimension L 400

A, while the (10,2), (10,3) and (11,2) SiCNTs
were placed in boxes withL 450

A. The boxes were pre-lled with water such that the density,
measured using the Voronoi tessellation method, was 1 g/cm
3
. Periodic boundary conditions were
imposed. After inserting the nanotubes in the middle of the simulation box, the energy of the
system was minimized, a step of the simulations that took about 5 ns. The Langevin thermostat
was used to increase the system's temperature until the desired T was reached, after which MD
simulation in the NPT ensemble was carried out to allow the water molecules to move in and
out of the nanotubes, which took about 15 ns. After equilibrium was reached, water molecules
outside the nanotubes were removed and the axis of the simulation box was aligned with that
of the nanotubes in order to simulate an innite nanotube. As the next step the system was
switched to the NP
zz
T ensemble with P
zz
1 atm, and MD simulations were performed for 30
ns or longer. Only the last 20 ns were used for collecting data and analyzing them.
Water molecules were represented by the TIP4P/Ice (Abascal et al., 2005) model, SiC was
modeled using the Terso potential (Terso, 1989), while the REBO (Brenner et al., 2002) po-
tential was utilized for the CNTs. The interactions between water and the nanotubes were rep-
resented by the Lennard-Jones (LJ) potential, whose parameters for SiCNTs were taken from
Malek and Sahimi (Malek and Sahimi, 2010), while the parameters of graphite were used for the
CNTs. The Lorentz-Berthelot mixing rules were employed to compute the LJ parameters of pairs
of atoms, while the partial charges for the SiC atoms were set to 0.6e, as reported based on ab
initio calculations (Mavrandonakis et al., 2003). The Coulombic interactions were computed with
the particle-particle-particle-mesh (Hockney and Eastwood, 1988) (PPPM) algorithm. To keep
the bonds and angles in the water molecules intact, the SHAKE algorithm (Ryckaert et al., 1977)
was used. The timestep was always 1 fs. All the simulations were carried out using the LAMMPS
66
(Plimpton, 1995) package, and the systems were visualized with the VMD (Humphrey et al., 1996)
package. To understand the eect of the partial charges on the HB formation-breaking process,
we also carried out simulations of water conned in CNTs at T 230 K with articially added
partial charges to the C atoms, with the values being 0.2e, 0.4e and 0.6e.
Several quantities were computed to address the question of whether the aforementioned rst-
order transition occurs. The dipole angle , and the axial and radial distribution functions were
computed in order to identify the possible melting points. The axial diusion coecients D
z
were
also computed, as were the cage correlation function, in order to check for the changes in the
cages surrounding the water molecules. These quantities and snapshots of the system were used
to determine if water was in the liquid or solid state. The range of temperature simulated was
from 140 K to 250 K, which was increased by 10 K steps.
4.2.2 Dipole Orientation
The dipole orientation angle , which corresponds to the projection of the dipole vector on
the plane formed by the unrolled nanotube and the z axis (Mei et al., 2015), was computed in
both SiCNTs and CNTs over the aforementioned range of temperatures, and its histograms were
constructed accordingly. In this way we rst ascertain that any dierence between the dipoles
was exclusively due the eect of the SiCNT and CNT connement, and not an artifact caused by
the choice of the water model.
4.2.3 The Axial and Radial Distribution Functions
To conrm the possible melting temperature identied in Figure 1 below, we computed gpzq,
the axial distribution function,
gpzq
1
N
N
¸
i1
N
¸
j1;ij
xpzz
ij
qy; (4.1)
67
where z
ij
is the distance between water molecules i and j along the z axis, and N is the total
number of the molecules. We used the oxygen atoms to determine z
ij
and compute gpzq. To
check further the evidence provided by the dipole orientation and the axial distribution function,
we also computed the radial distribution function gprq, based on the water oxygen atoms.
4.2.4 The Diusion Coecient
We computed the diusion coecient D
z
in the axial direction by two methods. In one we
computed the velocity autocorrelation function (VACF) of the oxygen atoms along the z axis,
from which the axial diusion coecient D
z
, given by the standard equation,
D
z

1
N
N
¸
i1
»
8
0
xv
zi
ptqv
zi
p0qydt; (4.2)
was computed, where v
zi
ptq is the velocity of oxygen i along the z-axis at time t, and xy im-
plies an ensemble average. In the second method D
z
was estimated through the mean-squared
displacements (MSDs):
D
z

1
2
lim
tÑ8
x|z
i
ptqz
i
p0q|
2
y
t
; (4.3)
where z
i
ptq is the axial position of oxygen i at time t. The velocity of the center of mass motion
was accounted for when computing the VACF. If diusion is Gaussian, then the MSDs would
increase linearly with the time, and Eq. (3) yields the correct diusion coecient. This is indeed
the case, as the MSDs do become a linear function of time after about 4-5 nanoseconds.
4.2.5 The Cage Correlation Function
Another quantity that was computed in order to conrm the diusion results was the cage
correlation function (Rabani et al., 1999) Cptq. If water is in a frozen state, then, one would
expect a cage formed by the neighboring molecules to not change over time. Thus, as the freezing
68
point is approached, the slow-down in the dynamics can be quantied by keeping track of the
molecular environment and how it changes. We used oxygen atoms as the basis for calculating
Cptq, because it is the closest to the center of mass of H
2
O.
To computeCptq we used a generalized neighbor list (Rabani et al., 1999) in which the locations
of the neighbors of every oxygen atom were recorded. Given an oxygen atom, a list of other atoms
in the rst solvation shell, which corresponds to the average distance between the peaks of the
axial distribution function, provides an ecient and accurate way of describing the immediate
neighborhood of the oxygen. If the location of any oxygen atom at time t remains unchanged
relative to time t 0, then, Cptq 1, and the list remains unchanged. If, however, a water
molecules diuses outside its surrounding cage, the lists at times t andt 0 will not be identical,
in which case the average Cptq 1 at time t. Thus, the cage correlation function Cptq was
computed as a function of time t by averaging over all the oxygen atoms in the nanotubes and
constructing the list for each oxygen atom and including all the neighbors of each individual atom
i in the SiCNT with radial distance r
list
,
L
i
ptqrfpr
ij
qs; j 1; 2; ;N (4.4)
where N is the total number of atoms in the system, and
fpr
ij
q
$
'
'
&
'
'
%
1 r
ij
¤r
list
0 otherwise
(4.5)
The computations take considerable time, if r
list
is large. In any case Cptq is given by
Cptq
xL
i
p0q L
i
ptqy
xL
2
i
p0qy
: (4.6)
69
To computeCptq, we built the neighbor lists using a cuto ofr
list
2:85

A, which is the distance
between the peaks in the gpzq.
4.2.6 Hydrogen Bond Correlation Function
A practical way to study the formation-breaking process of the HB networks is through com-
paring the intermittent HB correlation functions
8
(HBCF) in the nanotubes. To compute the
HBCF, a geometrical criterion was used (Barati Farimani and Aluru, 2016):
C
HB
ptq
xhptq hp0qy
xhp0q
2
y
; (4.7)
wherehptq represents the number of the HBs at timet, withhp0q being the corresponding number
at t 0. By denition, hp0q 1. Typically, the HBCF decays completely after only a few
picoseconds, much faster than the other correlation functions. It measures the probability of a
pair of water molecules being bonded at timet 0 and remaining bonded at timet, while ignoring
possible breaking within that time frame. In other words, the HBCF provides insight into how
persistent is a HB network. As stated earlier, we also computed the HBCF of water conned in
CNTs, by (articially) placing partial charges onto the C atoms whose values were 0.2e, 0.4e and
0.6e.
4.3 Results and discussion
4.3.1 Dipole orientation
The frequency distributions (or histogram) of the dipole orientations for water in the two
CNTs are shown in Figure 4.1. The most frequent dipole angle in both CNTs is about 30

(and
150

), which is close to those reported previously (Foroutan et al., 2020; Zeng et al., 2018). We
note that the references Foroutan et al. and Zeng et al. used the water models TIP3P and SPC/E,
70
respectively, which are not the same as what we used. This indicates that the two models and the
TIP4P/Ice that we use produce the same highest peak of the dipole angle . Moreover, having
computed the histograms in both CNTs and SiCNTs makes it possible to compare the dipole
angles of water in both nanotubes using the same water model.
0 20 40 60 80 100 120 140 160 180
T (K)
273
0 20 40 60 80 100 120 140 160 180
(degrees)
0.00
0.02
0.04
0.06
0.08
0.10
Abs. frequency
T (K)
273
(10,5) (11,3)
Figure 4.1: Histogram of the dipole angle of water conned in CNTs. The highest peaks indicate
the most frequent dipole angle, which is around 30

.
The corresponding histogram of dipole orientations for water in the SiCNTs are presented in
Figure 4.2. The most frequent dipole angle in all the nanotubes at all the temperatures is 40

(and
140

). This implies that for water in SiCNTs is independent of the diameter of the nanotubes
and the type of ice formed. The partial charges of Si and C, as well as the LJ parameters are
responsible for the orientation of the dipoles.
Figure 4.2 indicates a decrease in the height of peaks as the temperatureT increases. There is
a sharp decrease in the heights onceT is high enough, indicating that the ordered structure of the
water distribution in the nanotubes begins to disappear, hence signaling a phase transition, which
is usually interpreted as signaling a rst-order phase transition (Raju et al., 2018), to which we
shall return shortly. Nevertheless, after crossing the melting point, most of the water molecules
still favor the initial dipole angle angle of 40

. The MD simulation indicated (Cobe~ na Reyes et al.,
71
0.00
0.05
0.10
0.15
0.20
Abs. frequency
T (K)
140
160
180
190
T (K)
160
180
190
210
T (K)
150
170
180
200
0 20 40 60 80 100 120 140 160 180
0.00
0.05
0.10
0.15
0.20
T (K)
170
190
200
220
0 20 40 60 80 100 120 140 160 180
(degrees)
T (K)
170
190
200
220
0 20 40 60 80 100 120 140 160 180
T (K)
180
200
210
230
(9,2) (6,6) (8,4)
(10,2) (10,3) (11,2)
Figure 4.2: Histogram of the dipole angle of water conned in SiCNTs. The sharp drops in the
peaks indicate the possible melting points. The most frequent dipole angle is close to 40

2018) that this is due to the existence in the lower end of the temperature range of an icelike and
ordered, but dynamic, water. Therefore, although there is translational motion, there is still some
order in the distribution of the water molecules (see also below). Note that the relaxation times
in SiCNTs are up to two orders of magnitude larger than the bulk values (Cobe~ na-Reyes and
Sahimi, 2020a), which is consistent with the preferred position of the dipoles even after crossing
the melting point.
4.3.2 Ice types
Given the results so far, we proceeded to examine the type of ice-nanotubes formed. We
identied trigonal, square and pentagonal ice-nanotubes. That is, if we look into the nanotubes
through one end and along their axis, we see trigonal, square and pentagonal rings of water
molecules. Figure 4.3 presents snapshots of the shapes. The trigonal rings were formed in the
72
(9,2) nanotube, the square type was in the (6,6), (8,4) and (10,2) nanotubes, while the pentagonal
rings were formed in the (10,3) and (11,2) nanotubes. Chiashi et al. (Chiashi et al., 2019) found
that the same type of ice can form in various CNTs. This implies that water forms ice-nanotubes
of the same type, not only in CNTs, but also in SiCNTs, hence exhibiting a sort of universal
behavior.
(9,2) (6,6) (8,4)
(10,2) (10,3) (11,2)
Figure 4.3: The structure of the various types of ice in the SiCNTs.
4.3.3 The Radial and Axial Distribution Functions
Similar to the dipoles' frequency distributions shown in Figure 4.1, if there is melting, the
peaks in the axial distribution function (ADF) gpzq should also exhibit a sharp decrease, when
the temperature rises above the transition point. Figure 4.4, which presents the results for gpzq
in all the SiCNTs, does indicate this. The distance between two neighboring peaks is around
2:85

A, representing the distance between neighboring rings of water molecules that are formed
in the nanotubes (see also below), and is the same for all the nanotubes at all the temperatures.
73
g(z)OO
0
1
2
3
4
5
T (K)
140
160
180
190
T (K)
160
180
190
210
T (K)
150
170
180
200
0 2 4 6 8 10 12 14
0
1
2
3
4
5
T (K)
170
190
200
220
0 2 4 6 8 10 12 14
z(Å)
T (K)
170
190
200
220
0 2 4 6 8 10 12 14
T (K)
180
200
210
230
(9,2) (6,6) (8,4)
(10,2) (10,3) (11,2)
Figure 4.4: The axial distribution function of water in the SiCNTs. Similar to Fig. 4.1, a sudden
decrease in the peaks indicate a rst-order transition.
The results for the radial distribution function (RDF) gprq for all the nanotubes are shown in
Figure 4.5. They support the existence of a transition as the temperature rises. An interesting
nding is that some of the diagrams share some similarities. For example, the second peak in the
(6,6) nanotube is located at a distance shorter thanr 4

A, whereas the second peak in the (8,4)
nanotube is located at r 4

A. Both structures contain squared ice, but the smaller diameter of
the (6,6) nanotube has the water molecules slightly closer compared to the (8,4) nanotube. The
(10,2) nanotube also produces a squared ice, but due to its larger diameter the third peak develop
beyond r 5

A, which is due to a little more space being available among the molecules. The
same peak in the (6,6) and (8,4) nanotubes develops before r 5

A.
The (10,3) and (11,2) nanotubes also share a similar shape, as both produced a pentagonal
ice. There is, however, a small third peak in the (10,3) nanotube, which is too small in the (11,2)
nanotube. This may be because some of the rings in the pentagonal ice formed in the (10,3)
74
0
2
4
6
8
T (K)
160
180
190
210
0
1
2
3
4
5
T (K)
140
160
180
190
0
2
4
6
8
T (K)
150
170
180
200
0 1 2 3 4 5 6
0
2
4
6
8
T (K)
170
190
200
220
0 1 2 3 4 5 6
0
2
4
6
8
T (K)
180
200
210
230
g(r)OO
(9,2) (6,6)
(
(11,2) (10,2)
0 1 2 3 4 5 6
0
2
4
6
8
T (K)
170
190
200
220
r(Å)
(10,3)
(8,4)
Figure 4.5: The radial distribution function of water in the SiCNTs. Similar to Figs. 4.2 and 4.4,
a sudden decrease in the peaks indicate a rst-order transition.
nanotube surround a central water molecule, possibly due to the space restrictions. The smallest
nanotube, the (9,2), produces a unique shape due to its trigonal ice. Note that the second peak is
slightly bifurcated. That this happens is presumably due to dierences in the sizes of contiguous
trigonal rings, which are the result of vibrations in the nanotubes during the simulations.
4.3.4 The Diusion Coecient
After determining the possible temperatures where the ice-nanotubes melt, we computed the
diusion coecientD
z
. The results are presented in Figure 4.6, where we show the temperature-
dependence of the diusion coecient, computed by the VACF and the MSDs. The results
obtained by the two methods are similar, but the transition temperatures are clearer in the
results obtained through the VACFs. This is probably due to the fact that when computing D
z
using the MSDs, the rst several nanoseconds of about 20 ns were discarded in order to use only
75
the linear region of the dependence of the MSDs on time. The integration of the the entire VACF
does not suer from this drawback.
140 160 180 200 220 240 260
SiCNT
(9,2)
(6,6)
(8,4)
(10,2)
(10,3)
(11,2)
140 160 180 200 220 240 260
10
10
10
9
10
8
10
7
D
z
(cm
2
/s)
SiCNT
(9,2)
(6,6)
(8,4)
(10,2)
(10,3)
(11,2)
T (K)
Figure 4.6: Temperature-dependence of the axial diusion coecient D
z
, using the VACF (left)
and the MSDs (right). Arrows at lower T indicate the possible melting point, while those at
higher T indicate the fragile-to-strong transition.
As Figure 4.6 indicate, there is a clear transition in the diusion coecient with increasing
temperature. Water in the (9,2) SiCNT does not exhibit signicant changes in its diusivity
below T 180 K, but D
z
begins to increase steadily up to T 210 K, where a jump in the
curve develops. The rst change is related to the aforementioned phase transition. The sharp
and sudden upward change in the values ofD
z
is indicative of a rst-order transition. The second
change inD
z
pTq is presumably related to the fragile-to-strong (FTS) transition, because it occurs
close to the reported FTS transition temperature of T 230 K (Dubey et al., 2019; Shi et al.,
2018; Khademi and Sahimi, 2016; Mallamace et al., 2010). Similar behavior is seen in all the other
nanotubes. In the (6,6) and (8,4) nanotubes the rst-order transition occurs at around T 180
K, whereas the FTS transition temperature is at around T 210 K. The possible melting point
for the (10,2) tube is 180 K, and for the (10,3) is at T 190 K. The phase change in the (11,2)
occurs at T 210 K
76
4.3.5 The Cage Correlation Function
The computed Cptq for the six SiCNTs are showed in Figure 4.7. If there is diusion, one
expects the local environment around a water molecule to continuously change and, hence, Cptq
should decrease as the temperature increases. This occurred in all the nanotubes, but at dierent
temperatures, which is consistent with the results for the diusivity D
z
, the axial and radial
distribution functions, and the histogram of the dipole angle .
0.7
0.8
0.9
1.0
C(t)
T (K)
140
160
180
190
T (K)
160
180
190
210
T (K)
150
170
180
200
0.7
0.8
0.9
1.0
T (K)
170
190
200
220
time (ns)
T (K)
170
190
200
220
T (K)
180
200
210
230
(9,2) (6,6) (8,4)
(10,2) (11,2) (10,3)
10
6
10
5
10
4
10
3
10
2
10
1
10
0
10
1
10
6
10
5
10
4
10
3
10
2
10
1
10
0
10
1
10
6
10
5
10
4
10
3
10
2
10
1
10
0
10
1
Figure 4.7: Cage correlation function Cptq of the oxygen atoms. The decay of Cptq occurs at
temperatures that are consistent with those in Figs. 4.2, 4.4 and 4.5.
4.3.6 Strength of Hydrogen-Bond Networks
If the melting points in SiCNTs are lower than those in CNTs, it is likely that the HB network
in the former is weaker than in the latter, which may explain the depression in the freezing point
and the dierence between values of dipole orientation in the SiCNTs and CNTs. The HB
77
network is disturbed by the presence of partial charges in the nanotube, causing the ice structure
to be more unstable as the temperature increases. Indeed, it has been shown that the partial
charges do aect the average number of HBs per water molecules (Mei et al., 2015; Pokropivnyi
and Silenko, 2006).
We computed the HBCF in the SiCNTs at the melting temperature estimated earlier. The
results are presented in Figure 4.8. We also computed the same in CNTs at 273 K. The HBCF in
all the SiCNTs exhibit faster decay than those in the CNTs, even though the HBCF in the former
was computed at a much lower temperature than in the latter. This implies that even at lower
temperatures the HBCF is weaker in SiCNTs than in CNTs at higher temperatures. This implies
that it is likely that bifurcated HBs are formed (Foroutan et al., 2020), which is essentially a
step in the reorientation of the water molecules.
99
Thus, weak HB network is likely a cause of the
depression in the melting point in SiCNTs. We also note that in other types of nanotubes, such as
boron-nitrite nanotubes (Zeng et al., 2018), or when ionic species are present (Hernandez-Rojas
et al., 2012), the HBCF behaves in a completely dierent way than what we presented in Figure
4.8.
0.000 0.005 0.010 0.015 0.020 0.025 0.030
time (ns)
0.0
0.2
0.4
0.6
0.8
1.0
1.2
C
HB
(t)
(9,2)
(6,6)
(8,4)
(10,2)
(10,3)
(11,2)
(10,5)
(11,3)
Figure 4.8: Hydrogen-bond correlation function C
HB
of water in SiCNTs (circles) and CNTs
(squares).
78
4.3.7 Eect of Partial Charges.
It is important to understand how the HBCF changes when the partial charges of the C atoms
in the CNT are articially changed. As mentioned earlier, we carried out MD simulations in
which we changed the partial charges of the C atoms. The results, computed at T 273 K, are
presented in Figure 4.9. They indicate that the HBCF decays faster as partial charges increase,
i.e., the HB network does not last for too long compared to HBCF with the actual partial charges
that decays slower, implying that the melting point is lower than under the actual charges. We
note once again that even at higher temperatures, the decay of the HBCF with the altered partial
charges is slower, or at most similar to, than the correlation function in the SiCNTs. This is
because the water molecules interact more strongly with the charged C sites, which disrupt the
HB networks, re
ected in the resulting HBCF.
0.000 0.005 0.010 0.015 0.020 0.025 0.030
t (ns)
0.0
0.2
0.4
0.6
0.8
1.0
C
HB
(t)
Charge
0.0
0.2e
0.4e
0.6e
CNT
(10,5)
(11,3)
CNT
(10,5)
(11,3)
Figure 4.9: Eect of (articially) placed partial charges on the hydrogen bond correlation function
C
HB
of water in the CNTs atT 273 K. The rst two upper curves are the same as those in Fig.
4.8.
79
4.4 Conclusions
Extensive MD simulations were carried out in order to determine the melting point of water
conned in several SiCNTs of various sizes. The computed melting temperatures, which depend
on the size of the nanotubes, turn out to be signicantly lower than those reported by Chiashi
et al. (2019) for CNTs. The preferred dipole angle was computed to be 40

, which is also lower
than 30

for water in CNTs. Moreover, the HB network is weaker in SiCNTs than in CNTs,
which is presumably a factor that contributes to lowering the melting point in SiCNTs. The
computed diusivities exhibit not only the phase change, but also the fragile-to-strong transition
of supercooled water. The temperatures at which the FTS occurs varies slightly among the nan-
otubes. These results were supported by the computed axial and radial distribution functions, the
changes in the dipole angle distribution, and a comparison between the HB correlation functions
in SiCNTs and CNTs. We also studied the eect of the partial charges on the HBCF by changing
them in the CNTs, demonstrating that as the charges increase, the HBCF decays faster, hence
signaling a lower melting point. Distinct types of ice-nanotubes, including trigonal, square and
pentagonal, are also formed.
80
Chapter 5
Chemical potential and phase change of water in CNTs
using RexPoN potential
5.1 Introduction
Water under strong connement exhibits a very distinct behavior compared with the bulk
(Zhao et al., 2020; Yang and Guo, 2020; Cuadrado-Collados et al., 2020; Cobe~ na-Reyes and Sahimi,
2020b; Khademi and Sahimi, 2016). Its dynamic static and thermodynamic properties (Zhao et al.,
2020; Cobe~ na-Reyes and Sahimi, 2020a; Mondal and Bagchi, 2019; K ohler and Gavazzoni, 2019)
shift to values that can be orders of magnitude dierent from the value of bulk water (Cobe~ na-
Reyes and Sahimi, 2020b; Chiashi et al., 2019; Agrawal et al., 2017a). For instance, the freezing
point of conned water can drastically shift to lower or higher values (Cobe~ na-Reyes and Sahimi,
2020b; Chiashi et al., 2019; Agrawal et al., 2017a; Chakraborty et al., 2017; Khademi et al., 2015)
and properties such as viscosity under connement or supercooling can adopt a non-Newtonian
regime (Cobe~ na-Reyes and Sahimi, 2020a; de Almeida Ribeiro and de Koning, 2020).
Takaiwa et al. (2008) were the rst to report the phase diagram of water conned in CNTS
at atmospheric pressure. They used the TIP5P water model (Mahoney and Jorgensen, 2000) and
observed ice formation based on the change of enthalpy. Interestingly, in some of the nanotubes,
an abrupt phase change was observed while in others a gradual freezing occurred.
81
Later, Agrawal et al. (2017a) used Raman spectroscopy to nd that at a specic diameter of
nanotube water can exhibit an ice-like behavior at temperatures that were close to boiling point.
More recently, Chiashi et al. (2019) used photoluminescence spectroscopy to study the melting
point of water in CNTs. The results reported indicate that the change in the melting points was
not as drastic as previously thought.
Modeling water in CNTs is a dicult challenge. The models that are currently available are
designed to replicate some (but usually not all) properties of water. There are several popular
models used for this purpose such as the mW (Molinero and Moore, 2009), TIP4P-2005 (Abascal
and Vega, 2005) and the TIP4P/Ice (Abascal et al., 2005). The mW can reproduce accurately
the density at room temperature and the maximum density (Molinero and Moore, 2009) as well.
Similarly, it accurately reproduces the heat of vaporization and the formation of hexagonal ice
(Ih). The diusion coecient, however, is not one of the strengths of the model. The model
TIP4p-2005 accurately predicts the diusion coecient compared to other models such as TIP4P
and TIP5P at T=298 K. Additionally, the TIP4P-2005 can reproduce accurately the thermal
expansion coecient and the density of several ice polymorphs. The melting point is, however,
o by about 20

C degrees. The TIP4P/Ice was designed to study ice phases and to predict
accurately the bulk point and the formation of ice phases at dierent pressures. The TIP5P
predicts the melting point accurately too, but the polymorph produced is Ice II instead of Ih.
Nonetheless, reactive force elds have been recently used to study the properties of bulk and
conned water. Raju et al. (2018) used ReaxFF (Chenoweth et al., 2008) to study the phase
change of water in nanostructures. They lled up the nanotubes using the chemical potential
obtained from the particle insertion method proposed by Widom (1963) and showed that in some
cases the melting point was higher than 273 K. Naserifar and Goddard III (2019) recently used
the newly RexPoN (Naserifar and Goddard III, 2018) force eld to study properties of bulk water.
Their results indicate that the anomalies in water can be explained by a transformation of the
hydrogen bond network from 1D to 2D. This last force eld, the RexPoN potential, was formulated
82
using only quantum mechanics (QM) data and with no empiric formulae. The RexPoN model is
capable to reproduce almost all of the water properties in a wide range of temperatures including
the melting point T
m
= 273.3 K (expr = 273.15 K) and properties at 298 K: Hvap = 10.36
kcal/mol (expr = 10.52), density = 0.9965 gr/cm
3
(expr = 0.9965), entropy = 68.4 (J/mol)/K
(expr = 69.9), dielectric constant = 76.1 (expr = 78.4). And more recently, RexPoN was updated
with corrections for other elements such as H, C, N, F, Cl, Br, I, P and the noble gases (Naserifar
et al., 2019). RexPoN, however, has not been used to model a system such as water + CNT. In
this chapter we aim to study the phase change of water conned in CNTs using this model. In
this chapter we attempt to address the question of what is the phase diagram of water conned
in CNTs using a reactive force eld as RexPoN.
First we started computing the chemical potential of water in a range of temperatures using
the Widom particle insertion method (Widom, 1963). Then we proceeded with MC simulations
to obtain structures that only depend on the chemical potential. Finally we used MD simulations
to compute several properties that would help us to see the phase diagram of water.
5.2 Computational Methods
5.2.1 Description of the RexPoN potential
The RexPoN potential is composed by several energy terms. First, the long range interactions
such as the van der Waals forces, which comprises the London dispersion and the Pauli repulsion
are derived using QM. Second, the electrostatics are modeled using the PQEq (Naserifar et al.,
2017) method which reproduces accurately QM polarization energies of small molecules. Thus,
RexPoN supports charge polarization. Additionally, the model accounts for the hydrogen bond
(HB) corrections and is capable to determine bond and angle valence terms by using the Coupled
83
Cluster Single Double Triple CCSD(T) (Shank et al., 2009) on the monomer and dimer of water.
Based on this The total energy of the system is expressed as:
E
total
E
nonbond
E
valence
(5.1)
whereE
nonbond
comprises the non-bond and bonded energy terms respectively. Then, these 2
terms are dened by:
E
nonbond
E
elect
E
vdW
E
HB
(5.2)
E
valence
E
bond
E
angle
(5.3)
The van der Waals energy E
elect
was obtained using the equations of state of solid phases of
oxygen and hydrogen based on density functional theory while the electrostatic energy E
elect
is
described using the PQEq formulation as mentioned before.
The covalent bond energy E
bond
used the bond order (BO) formulation which uses the same
sigmoid function as in ReaxFF (Van Duin et al., 2001). The angle energy E
angle
is dened as
a cosine angle function that was tted to the potential energy surface of the water monomer as
expressed by Shank et al. (2009).
To compute the HB energy E
HB
the geometric orientation of the donor and acceptor water
molecules must be accounted. Hence, the interaction of a lone-pair and of several training sets
were used to account for possible conguration regarding the orientation of possible dimers.
5.2.2 Chemical potential of the RexPoN water model
To compute the chemical potential the particle-insertion method proposed by Widom Widom
(1963) whose algorithm in LAMMPS is explained in Frenkel and Smit (2001) was used. In the
84
Widom method, the excess chemical potential is computed using the dierence between the total
energy of the system and the interaction energy of the system with the new added particle. This
is mathematically expressed as Frenkel and Smit (2001):

ex
k
B
T ln
»
ds
N1
xexprUsy
N
; (5.4)
where
ex
k
B
, s, and U are the excess chemical potential, the Boltzman constant, the
spatial coordinates, the reciprocal temperature (1{k
B
T ) and the dierence of potential energy
between the systems with N 1 and N particles respectively.
The above equation comes from the denition of the Helmholtz energy and the denition of
the chemical potential:

id

ex
(5.5)
where the subscriptsid andex stand for ideal gas and excess contributions respectively. Equation
5.2 can be rewritten as:
k
B
T lnpQ
N1
{Q
N
q; (5.6)
where Q is the partition function. The complete expression for the chemical potential using
the partition function is:
k
B
T ln

V{
d
N 1

k
B
T ln
#
³
ds
N1
exprUps
N1
qs
³
ds
N
exprUps
N
qs
+
; (5.7)
where V , d and , are the volume, the dimensionality and the thermal de Broigle wavelength
respectively (Frenkel and Smit, 2001).
It must be noted that the \inserted" particles are never accepted. In other words, the molecules
that are added to compute the interaction energy of the N 1 particle are deleted after the
85
computation. One can think of these molecules as \virtual". This is done to ensure that the total
number of molecules is not altered during the simulation.
To compute the
ex
of pure water a box of L 15

A containing 256 water molecules was
used. The timestep was 1 fs. All water molecules were modeled using the RexPoN potential.
The range of temperatures for all simulations was 140K - 340K. We used the NPT ensemble
and set P = 1 atm. After a minimization and an equilibration step of about 2 ns to reach the
desired temperature, the particles were \inserted" at a rate of 2 particles every 1000 timesteps.
This procedure was performed for 10 ns allowing enough time for sampling dierent insertion
locations.
5.2.3 Filling up the CNTs with water
In the previous chapters, the method used to ll up the nanotubes was to \soak" the tubes
in a box full of water. While this method is totally valid, it can be inecient, specially because
we are not interested in the water molecules located outside the nanotube. Furthermore, the
RexPoN potential is expensive because of its reactive nature and simulating such a large system
will require additional computational time and resources. As mentioned in Chapter 1, a more
ecient approach is to use the Grand Canonical Monte Carlo (GCMC) ensemble to ll up the
nanotubes. In this Chapter we used three nanotubes: (11,2), (12,5) and (20,6) with diameters of
9:56, 11:53 and 18:59

Arespectively. The lengths of the tubes were: 416, 258 and 152

A. The valid
MC moves were insertion, deletion and rotation of the water molecules added in to the nanotube.
The MC moves were performed until a constant number of water molecules were present inside
the nanotube. Basically, the nanotube gets water molecules added, removed or rotated until the
chemical potential of the water molecules inside the nanotube is equal to the chemical potential
of the ctitious gas reservoir.
Alternate processes were used: MC moves followed by time integration in theNP
zz
T ensemble
with P
zz
1 atm. This is necessary to allow for reorganization of the molecules inside the
86
nanotube. At least 20 cycles were performed. Each MC process comprises 1000 MC moves. Then
the time integration for 1000000 steps was performed. The time step was always 1.0 fs.
5.2.4 Molecular dynamics of water in CNTs using RexPoN
Once obtained the structures of water + CNT, MD simulations were extensively used to
observe the change in the potential energy (PE) and the radial distribution function (gprq) of the
water molecules. TheNP
zz
T ensemble was used such thatP
zz
was set to 1 atm as in the previous
step. The dynamics were run for 10 ns using a timestep of 1 fs. A rst equilibration step of 4 ns
was used before the run of 10 ns. The data of the (PE) and the gprq were collected at intervals of
0.01 ns.
5.3 Results and discussion
The data corresponding to the excess chemical potential
ex
of water and supercooled water
is shown in Figure 5-1. The values of
ex
increases as the Temperature increase. This is obvious,
since
ex
is basically the contribution of the non-ideal behavior of water. The values of
ex
were used to ll up the nanotubes at each temperatures to generate the structures. For practical
purposes we used supercooled water instead of ice. This is important because the ice formed inside
the nanotubes is not the same as the bulk ice phases known. Finally, the values of liquid water of

ex
are in good agreement with theoretical values reported by others (D om ot or and Hentschke,
2004; Rahbari et al., 2018).
After the MC simulation 322, 408 and 507 water molecules at T = 300 K were found inside
the nanotubes. Using Voronoi tessellation these numbers are equivalent to roughly 0:97 0:98
g/cm
3
. These results are similar to the density found in SiCNTs as shown in Chapter 4.
Inside the nanotubes three dierent type of ices were found: Trigonal, hexagonal and octagonal
ices in the (11,2), (12,5) and (20,6) CNTs respectively. Similarly, Takaiwa et al. (2008) found
87
-32
-30
-28
-26
-24
-22
100 150 200 250 300 350
μ
ex
(kJ/mol)
T (K)
Figure 5.1: Excess chemical potential of water
(11,2) (12,5)
(20,6)
Figure 5.2: Ice phases observed in the three CNTs.
dierent ice phases. However, they used the TIP5P water model which produces Ice II instead of
Ih (Abascal and Vega, 2005). RexPoN produces hexagonal ice instead.
88
-45
-40
-35
-30
-25
100 150 200 250 300 350
(11,2)
(12,5)
(20,6)
Potential Energy (kJ/mol)
T (K)
Figure 5.3: Potential energy of water in the three CNTs
0
1
2
3
4
5
6
7
8
0 1 2 3 4 5
200
270
g(r)
0
2
4
6
8
10
0 1 2 3 4 5 6
200
290
r(Å)
0
1
2
3
4
5
6
7
0 1 2 3 4 5 6 7
200
330
(11,2) (12,5)
(20,6)
Figure 5.4: Radial distribution function gprq of water below and above melting point.
89
The ice phases at T = 260K are showed in Figure 5-2. In Chapter 4 we reported the presence
of similar ices. Nonetheless, the melting points in SiCNTs are very dierent than in CNTs. The
change in the PE as a function of temperature is presented in Figure 5-3. The abrupt change in
PE is about 6 kJ/mol which is close to the value of the heat of fusion of water. Thus, RexPoN
possibly also reproduces correctly the experimental value of the heat of fusion.
The melting points in each of the nanotubes are 250, 270 and 310 K for (11,2), (12,6) and
(20,6) respectively. These values are not as high as reported by Agrawal et al. (2017a). Instead,
they are more similar to those reported by Chiashi et al. (2019). It is worth noting that the
results in this Chapter, in contrast to the previous chapters, a single force eld that models the
whole system was used. This is an advantage that can lead to more accurate results. It would
be interesting to model a system of SiCNT + water using RexPoN and to compare the values
obtained in Chapter 4.
Information about the ice phases and the melting point of water in the CNTs can also be
obtained by computing the radial distribution function. The results about are shown in Figure 5-
4. There is a clear dierence between the the g(r) of both temperatures. In all three temperatures,
the drop in the height of the peaks is consistent with the values obtained in the potential energy
diagram.
5.4 Conclusion
In this Chapter the excess chemical potential of the RexPoN water model was computed. We
hope that these result will be useful to others when modeling adsorption of water in nanostructures.
Moreover, the melting point in three dierent CNTs does not directly depend on the diameter
and that the change in the melting point is mild compared to previous results. Finally, dierent
ice phases were found with dierent melting points. This is in agreement with the results shown
in previous Chapters.
90
Bibliography
Abascal, J. L. and Vega, C. (2005). A general purpose model for the condensed phases of water:
Tip4p/2005. The Journal of chemical physics, 123(23):234505.
Abascal, J. L. F., Sanz, E., Fernandez, R. G., and Vega, C. (2005). A potential model for the
study of ices and amorphous water: Tip4p/ice. Journal of Chemical Physics, 122(23):9.
Abbas, A., Al-Amer, A. M., Laoui, T., Al-Marri, M. J., Nasser, M. S., Khraisheh, M., Atieh,
M. A., et al. (2016). Heavy metal removal from aqueous solution by advanced carbon nanotubes:
critical review of adsorption applications. Separation and Purication Technology, 157:141{161.
Abdelhalim, A., Winkler, M., Loghin, F., Zeiser, C., Lugli, P., and Abdellah, A. (2015). Highly
sensitive and selective carbon nanotube-based gas sensor arrays functionalized with dierent
metallic nanoparticles. Sensors and Actuators B: Chemical, 220:1288{1296.
Abdullah, M. and Kamarudin, S. K. (2017). Titanium dioxide nanotubes (tnt) in energy and
environmental applications: An overview. Renewable and Sustainable Energy Reviews, 76:212{
225.
Agarwal, A., Bakshi, S. R., and Lahiri, D. (2016). Carbon nanotubes: reinforced metal matrix
composites. CRC press.
Agrawal, K. V., Shimizu, S., Drahushuk, L. W., Kilcoyne, D., and Strano, M. S. (2017a). Ob-
servation of extreme phase transition temperatures of water conned inside isolated carbon
nanotubes. Nature Nanotechnology, 12(3):267.
Agrawal, K. V., Shimizu, S., Drahushuk, L. W., Kilcoyne, D., and Strano, M. S. (2017b). Ob-
servation of extreme phase transition temperatures of water conned inside isolated carbon
nanotubes. Nature nanotechnology, 12(3):267.
Alba-Simionesco, C., Coasne, B., Dosseh, G., Dudziak, G., Gubbins, K. E., Radhakrishnan, R.,
and Sliwinska-Bartkowiak, M. (2006). Eects of connement on freezing and melting. Journal
of Physics-Condensed Matter, 18(6):R15{R68.
Angell, C. A., Bressel, R. D., Hemmati, M., Sare, E. J., and Tucker, J. C. (2000). Water and its
anomalies in perspective: tetrahedral liquids with and without liquid-liquid phase transitions.
Physical Chemistry Chemical Physics, 2(8):1559{1566.
Ansari, N., Onat, B., Sosso, G. C., and Hassanali, A. (2020). Insights into the emerging networks
of voids in simulated supercooled water. The Journal of Physical Chemistry B, 124(11):2180{
2190.
Azamat, J., Khataee, A., and Joo, S. W. (2014). Separation of a heavy metal from water through
a membrane containing boron nitride nanotubes: molecular dynamics simulations. Journal of
Molecular Modeling, 20(10).
91
Bai, D. (2011). Size, morphology and temperature dependence of the thermal conductivity of
single-walled silicon carbide nanotubes. Fullerenes Nanotubes and Carbon Nanostructures,
19(4):271{288.
Banys, J., Kinka, M., Macutkevic, J., Volkel, G., Bohlman, W., Umamaheswari, V., Hartmann,
M., and Poppl, A. (2006). Eect of connement on the freezing-melting dynamics of water,
volume 514-516 of Materials Science Forum, pages 1255{1259.
Barati Farimani, A. and Aluru, N. R. (2016). Existence of multiple phases of water at nanotube
interfaces. The Journal of Physical Chemistry C, 120(41):23763{23771.
Barghi, S. H., Tsotsis, T. T., and Sahimi, M. (2014a). Hydrogen sorption hysteresis and superior
storage capacity of silicon-carbide nanotubes over their carbon counterparts. International
Journal of Hydrogen Energy, 39(36):21107{21115.
Barghi, S. H., Tsotsis, T. T., and Sahimi, M. (2014b). Hydrogen sorption hysteresis and superior
storage capacity of silicon-carbide nanotubes over their carbon counterparts. international
journal of hydrogen energy, 39(36):21107{21115.
Barghi, S. H., Tsotsis, T. T., and Sahimi, M. (2016a). Experimental investigation of hydro-
gen adsorption in doped silicon-carbide nanotubes. International Journal of Hydrogen Energy,
41(1):369{374.
Barghi, S. H., Tsotsis, T. T., and Sahimi, M. (2016b). Experimental investigation of hydro-
gen adsorption in doped silicon-carbide nanotubes. international journal of hydrogen energy,
41(1):369{374.
Berber, S., Kwon, Y.-K., and Tom anek, D. (2000). Unusually high thermal conductivity of carbon
nanotubes. Physical review letters, 84(20):4613.
Bernal, J. D. and Fowler, R. H. (1933). A theory of water and ionic solution, with particular
reference to hydrogen and hydroxyl ions. The Journal of Chemical Physics, 1(8):515{548.
Bertrand, C. E., Liu, K. H., Mamontov, E., and Chen, S. H. (2013a). Hydration-dependent
dynamics of deeply cooled water under strong connement. Physical Review E, 87(4).
Bertrand, C. E., Zhang, Y., and Chen, S. H. (2013b). Deeply-cooled water under strong conne-
ment: neutron scattering investigations and the liquid-liquid critical point hypothesis. Physical
Chemistry Chemical Physics, 15(3):721{745.
Bhadra, M., Roy, S., and Mitra, S. (2016). A bilayered structure comprised of functionalized car-
bon nanotubes for desalination by membrane distillation. ACS applied materials & interfaces,
8(30):19507{19513.
Bostick, D. and Berkowitz, M. L. (2003). The implementation of slab geometry for membrane-
channel molecular dynamics simulations. Biophysical journal, 85(1):97{107.
Brenner, D. W. (1990). Empirical potential for hydrocarbons for use in simulating the chemical
vapor deposition of diamond lms. Physical review B, 42(15):9458.
Brenner, D. W., Shenderova, O. A., Harrison, J. A., Stuart, S. J., Ni, B., and Sinnott, S. B.
(2002). A second-generation reactive empirical bond order (rebo) potential energy expression
for hydrocarbons. Journal of Physics-Condensed Matter, 14(4):783{802.
Briganti, G., Rogati, G., Parmentier, A., Maccarini, M., and De Luca, F. (2017). Neutron
scattering observation of quasi-free rotations of water conned in carbon nanotubes. Scientic
reports, 7(1):1{10.
92
Brovchenko, I., Geiger, A., and Oleinikova, A. (2003). Multiple liquid-liquid transitions in super-
cooled water. Journal of Chemical Physics, 118(21):9473{9476.
Budyka, M., Zyubina, T., Ryabenko, A., Lin, S., and Mebel, A. (2005). Bond lengths and
diameters of armchair single wall carbon nanotubes. Chemical Physics Letters, 407(4):266 {
271.
Cao, W., Huang, L., Ma, M., Lu, L., and Lu, X. (2018). Water in narrow carbon nanotubes:
roughness promoted diusion transition. The Journal of Physical Chemistry C, 122(33):19124{
19132.
Cao, Z., Peng, Y., Yan, T., Li, S., Li, A., and Voth, G. A. (2010). Mechanism of fast proton
transport along one-dimensional water chains conned in carbon nanotubes. Journal of the
American Chemical Society, 132(33):11395{11397.
Chakraborty, S., Kumar, H., Dasgupta, C., and Maiti, P. K. (2017). Conned water: structure,
dynamics, and thermodynamics. Accounts of chemical research, 50(9):2139{2146.
Chen, S. H., Gallo, P., Sciortino, F., and Tartaglia, P. (1997). Molecular-dynamics study of
incoherent quasielastic neutron-scattering spectra of supercooled water. Physical Review E,
56(4):4231{4243.
Chen, S.-H., Mallamace, F., Mou, C.-Y., Broccio, M., Corsaro, C., Faraone, A., and Liu, L.
(2006). The violation of the stokes{einstein relation in supercooled water. Proceedings of the
National Academy of Sciences, 103(35):12974{12978.
Chen, X., Cao, G., Han, A., Punyamurtula, V. K., Liu, L., Culligan, P. J., Kim, T., and Qiao, Y.
(2008). Nanoscale
uid transport: size and rate eects. Nano letters, 8(9):2988{2992.
Chenoweth, K., Van Duin, A. C., and Goddard, W. A. (2008). Reax reactive force eld for
molecular dynamics simulations of hydrocarbon oxidation. The Journal of Physical Chemistry
A, 112(5):1040{1053.
Chiashi, S., Saito, Y., Kato, T., Konabe, S., Okada, S., Yamamoto, T., and Homma, Y. (2019).
Connement eect of sub-nanometer dierence on melting point of ice-nanotubes measured by
photoluminescence spectroscopy. ACS nano, 13(2):1177{1182.
Chong, S.-H. and Ham, S. (2016). Anomalous dynamics of water conned in protein{protein and
protein{dna interfaces. The journal of physical chemistry letters, 7(19):3967{3972.
Chu, X.-Q., Kolesnikov, A., Moravsky, A., Garcia-Sakai, V., and Chen, S.-H. (2007). Observation
of a dynamic crossover in water conned in double-wall carbon nanotubes. Physical Review E,
76(2):021505.
Ciora, R. J., Fayyaz, B., Liu, P. K., Suwanmethanond, V., Mallada, R., Sahimi, M., and Tsotsis,
T. T. (2004). Preparation and reactive applications of nanoporous silicon carbide membranes.
Chemical Engineering Science, 59(22-23):4957{4965.
Cobe~ na Reyes, J., Kalia, R. K., and Sahimi, M. (2018). Complex behavior of ordered and icelike
water in carbon nanotubes near its bulk boiling point. The journal of physical chemistry letters,
9(16):4746{4752.
Cobe~ na-Reyes, J. and Sahimi, M. (2020a). Rheology of water in small nanotubes. Physical Review
E, 102(2):023106.
Cobe~ na-Reyes, J. and Sahimi, M. (2020b). Universal intrinsic dynamics and freezing of water in
small nanotubes. The Journal of Physical Chemistry C.
93
Coleman, J. N., Khan, U., Blau, W. J., and Gun'ko, Y. K. (2006). Small but strong: a review of
the mechanical properties of carbon nanotube{polymer composites. Carbon, 44(9):1624{1652.
Cuadrado-Collados, C., Majid, A. A., Martinez-Escandell, M., Daemen, L. L., Missyul, A., Koh,
C., and Silvestre-Albero, J. (2020). Freezing/melting of water in the conned nanospace of
carbon materials: Eect of an external stimulus. Carbon, 158:346{355.
Cummings, P. and Stell, G. (1981). Mean spherical approximation for a model liquid metal
potential. Molecular Physics, 43(6):1267{1291.
Cupane, A., Fomina, M., Piazza, I., Peters, J., and Schiro, G. (2014). Experimental evidence for
a liquid-liquid crossover in deeply cooled conned water. Physical Review Letters, 113(21).
D'Alessandro, A., Rallini, M., Ubertini, F., Materazzi, A. L., and Kenny, J. M. (2016). Inves-
tigations on scalable fabrication procedures for self-sensing carbon nanotube cement-matrix
composites for shm applications. Cement and Concrete Composites, 65:200{213.
Das, R., Ali, M. E., Abd Hamid, S. B., Ramakrishna, S., and Chowdhury, Z. Z. (2014). Carbon
nanotube membranes for water purication: A bright future in water desalination. Desalination,
336:97{109.
de Almeida Ribeiro, I. and de Koning, M. (2020). Non-newtonian
ow eects in supercooled
water. Physical Review Research, 2(2):022004.
de la Torre, J., Duque-Zumajo, D., Camargo, D., and Espa~ nol, P. (2019). Microscopic slip
boundary conditions in unsteady
uid
ows. Physical Review Letters, 123(26):264501.
Debenedetti, P. G. (1996). Metastable liquids: concepts and principles. Princeton University
Press.
Debenedetti, P. G. (2003). Supercooled and glassy water. Journal of Physics-Condensed Matter,
15(45):R1669{R1726.
Debenedetti, P. G. and Stillinger, F. H. (2001). Supercooled liquids and the glass transition.
Nature, 410(6825):259{267.
Dehaoui, A., Issenmann, B., and Caupin, F. (2015). Viscosity of deeply supercooled water
and its coupling to molecular diusion. Proceedings of the National Academy of Sciences,
112(39):12020{12025.
Dellago, C., Naor, M. M., and Hummer, G. (2003). Proton transport through water-lled carbon
nanotubes. Physical review letters, 90(10):105902.
Deserno, M. and Holm, C. (1998a). How to mesh up ewald sums. i. a theoretical and numerical
comparison of various particle mesh routines. The Journal of chemical physics, 109(18):7678{
7693.
Deserno, M. and Holm, C. (1998b). How to mesh up ewald sums. ii. an accurate error estimate for
the particle{particle{particle-mesh algorithm. The Journal of chemical physics, 109(18):7694{
7701.
Diallo, S. O. (2015). Pore-size dependence and characteristics of water diusion in slitlike micro-
pores. Physical Review E, 92(1):012312.
Diallo, S. O., Vlcek, L., Mamontov, E., Keum, J. K., Chen, J. H., Hayes, J. S., and Chialvo,
A. A. (2015). Translational diusion of water inside hydrophobic carbon micropores studied by
neutron spectroscopy and molecular dynamics simulation. Physical Review E, 91(2).
94
Dokmaj, T., Ibrahim, T., Khamis, M., Abouleish, M., and Alam, I. (2020). Chemically modied
nanoparticles usage for removal of chromium from sewer water. Environmental Nanotechnology,
Monitoring & Management, page 100319.
D om ot or, G. and Hentschke, R. (2004). Atomistically modeling the chemical potential of small
molecules in dense systems. The Journal of Physical Chemistry B, 108(7):2413{2417.
Donati, C., Glotzer, S. C., Poole, P. H., Kob, W., and Plimpton, S. J. (1999). Spatial correla-
tions of mobility and immobility in a glass-forming lennard-jones liquid. Physical Review E,
60(3):3107.
Dong, Y., Ma, L., Tang, C. Y., Yang, F., Quan, X., Jassby, D., Zaworotko, M. J., and Guiver,
M. D. (2018). Stable superhydrophobic ceramic-based carbon nanotube composite desalination
membranes. Nano letters, 18(9):5514{5521.
Dubey, V., Erimban, S., Indra, S., and Daschakraborty, S. (2019). Understanding the origin of
the breakdown of the stokes{einstein relation in supercooled water at dierent temperature{
pressure conditions. The Journal of Physical Chemistry B, 123(47):10089{10099.
Eatemadi, A., Daraee, H., Karimkhanloo, H., Kouhi, M., Zarghami, N., Akbarzadeh, A., Abasi,
M., Hanifehpour, Y., and Joo, S. W. (2014). Carbon nanotubes: properties, synthesis, puri-
cation, and medical applications. Nanoscale research letters, 9(1):393.
Ebrahimi, F., Maktabdaran, G. R., and Sahimi, M. (2020). Formation of a stable bridge between
two disjoint nanotubes with single-le chains of water. The Journal of Physical Chemistry B,
124(38):8340{8346.
Ebrahimi, F. and Pishevar, A. (2015). Dependence of the dynamics of spontaneous imbibition into
carbon nanotubes on the strength of molecular interactions. The Journal of Physical Chemistry
C, 119(51):28389{28395.
Ebrahimi, F., Ramazani, F., and Sahimi, M. (2018). Nanojunction eects on water
ow in carbon
nanotubes. Scientic reports, 8(1):1{10.
Elyassi, B., Sahimi, M., and Tsotsis, T. T. (2007). Silicon carbide membranes for gas separation
applications. Journal of Membrane Science, 288(1-2):290{297.
Elyassi, B., Sahimi, M., and Tsotsis, T. T. (2008). A novel sacricial interlayer-based method for
the preparation of silicon carbide membranes. Journal of Membrane Science, 316(1-2):73{79.
Eskhan, A., Banat, F., Abu Haija, M., and Al-Asheh, S. (2019). Synthesis of meso-
porous/macroporous microparticles using three-dimensional assembly of chitosan-functionalized
halloysite nanotubes and their performance in the adsorptive removal of oil droplets from water.
Langmuir, 35(6):2343{2357.
Evans, D. J. and Morriss, G. (1984). Nonlinear-response theory for steady planar couette
ow.
Physical Review A, 30(3):1528.
Fanetti, S., Lapini, A., Pagliai, M., Citroni, M., Di Donato, M., Scandolo, S., Righini, R., and
Bini, R. (2013). Structure and dynamics of low-density and high-density liquid water at high
pressure. The journal of physical chemistry letters, 5(1):235{240.
Farmahini, A. H., Shahtalebi, A., Jobic, H., and Bhatia, S. K. (2014). In
uence of structural
heterogeneity on diusion of ch4 and co2 in silicon carbide-derived nanoporous carbon. The
Journal of Physical Chemistry C, 118(22):11784{11798.
95
Fernandez, R. G., Abascal, J. L. F., and Vega, C. (2006). The melting point of ice i-h for
common water models calculated from direct coexistence of the solid-liquid interface. Journal
of Chemical Physics, 124(14).
Foroutan, M., Fatemi, S. M., and Shokouh, F. (2016). Graphene connement eects on melt-
ing/freezing point and structure and dynamics behavior of water. Journal of Molecular Graphics
& Modelling, 66:85{90.
Foroutan, M., Naeini, V. F., and Ebrahimi, M. (2020). Enhanced wettability of long narrow
carbon nanotubes in a double-walled hetero-structure: unraveling the eects of a boron nitride
nanotube as the exterior. Physical Chemistry Chemical Physics, 22(1):391{401.
Frenkel, D. and Smit, B. (2001). Understanding molecular simulation: from algorithms to appli-
cations, volume 1. Elsevier.
Furukawa, A. (2017). Onset of shear thinning in glassy liquids: Shear-induced small reduction of
eective density. Physical Review E, 95(1):012613.
Gallo, P., Rovere, M., and Chen, S. H. (2010). Dynamic crossover in supercooled conned water:
Understanding bulk properties through connement. Journal of Physical Chemistry Letters,
1(4):729{733.
Gallo, P., Sciortino, F., Tartaglia, P., and Chen, S. H. (1996). Slow dynamics of water molecules
in supercooled states. Physical Review Letters, 76(15):2730{2733.
Gao, H., Shi, Q., Rao, D., Zhang, Y., Su, J., Liu, Y., Wang, Y., Deng, K., and Lu, R. (2017).
Rational design and strain engineering of nanoporous boron nitride nanosheet membranes for
water desalination. The Journal of Physical Chemistry C, 121(40):22105{22113.
Ge, M., Li, Q., Cao, C., Huang, J., Li, S., Zhang, S., Chen, Z., Zhang, K., Al-Deyab, S. S.,
and Lai, Y. (2017). One-dimensional tio2 nanotube photocatalysts for solar water splitting.
Advanced Science, 4(1):1600152.
Gohardani, O., Elola, M. C., and Elizetxea, C. (2014). Potential and prospective implementation
of carbon nanotubes on next generation aircraft and space vehicles: A review of current and
expected applications in aerospace sciences. Progress in Aerospace Sciences, 70:42{68.
Gu, W. and Yushin, G. (2014). Review of nanostructured carbon materials for electrochemical
capacitor applications: advantages and limitations of activated carbon, carbide-derived carbon,
zeolite-templated carbon, carbon aerogels, carbon nanotubes, onion-like carbon, and graphene.
Wiley Interdisciplinary Reviews: Energy and Environment, 3(5):424{473.
He, Z., Zhou, J., Lu, X., and Corry, B. (2013). Ice-like water structure in carbon nanotube (8, 8)
induces cationic hydration enhancement. The Journal of Physical Chemistry C, 117(21):11412{
11420.
Hernandez-Rojas, J., Calvo, F., Bret on, J., and Gomez Llorente, J. (2012). Connement eects on
water clusters inside carbon nanotubes. The Journal of Physical Chemistry C, 116(32):17019{
17028.
Hockney, R. W. and Eastwood, J. W. (1988). Computer simulation using particles. crc Press.
Holten, V., Bertrand, C. E., Anisimov, M. A., and Sengers, J. V. (2012). Thermodynamics of
supercooled water. Journal of Chemical Physics, 136(9):18.
96
Horn, H. W., Swope, W., Pitera, J., Madura, J. D., Dick, T. J., Hura, G. L. B., and Head-Gordon,
T. (2004). Development of an improved four-site water model for bio-molecular simulations:
Tip4p-ew. Abstracts of Papers of the American Chemical Society, 228:U530{U531.
Horn, H. W., Swope, W. C., and Pitera, J. W. (2005). Characterization of the tip4p-ew water
model: Vapor pressure and boiling point. Journal of Chemical Physics, 123(19).
Hub, J. S., de Groot, B. L., Grubmuller, H., and Groenhof, G. (2014). Quantifying artifacts in
ewald simulations of inhomogeneous systems with a net charge. Journal of chemical theory and
computation, 10(1):381{390.
Hughes, Z. E., Shearer, C. J., Shapter, J., and Gale, J. D. (2012). Simulation of water trans-
port through functionalized single-walled carbon nanotubes (swcnts). The Journal of Physical
Chemistry C, 116(47):24943{24953.
Hummer, G., Rasaiah, J. C., and Noworyta, J. P. (2001). Water conduction through the hy-
drophobic channel of a carbon nanotube. Nature, 414(6860):188.
Humphrey, W., Dalke, A., and Schulten, K. (1996). Vmd: Visual molecular dynamics. Journal
of Molecular Graphics & Modelling, 14(1):33{38.
Iijima, S. (1991). Helical microtubules of graphitic carbon. Nature, 354(6348):56{58.
Ingebrigtsen, T. S. and Tanaka, H. (2018). Structural predictor for nonlinear sheared dynamics
in simple glass-forming liquids. Proceedings of the National Academy of Sciences, 115(1):87{92.
Jabbarzadeh, A., Harrowell, P., and Tanner, R. (2006). Crystal bridge formation marks the
transition to rigidity in a thin lubrication lm. Physical review letters, 96(20):206102.
Jin, J., Fu, L., Yang, H., and Ouyang, J. (2015). Carbon hybridized halloysite nanotubes for
high-performance hydrogen storage capacities. Scientic reports, 5:12429.
Jorgensen, W. L., Chandrasekhar, J., Madura, J. D., Impey, R. W., and Klein, M. L. (1983).
Comparison of simple potential functions for simulating liquid water. The Journal of chemical
physics, 79(2):926{935.
Kafshgari, M. H., Mazare, A., Distaso, M., Goldmann, W. H., Peukert, W., Fabry, B., and
Schmuki, P. (2019). Intracellular drug delivery with anodic titanium dioxide nanotubes and
nanocylinders. ACS applied materials & interfaces, 11(16):14980{14985.
Karimi, M., Solati, N., Ghasemi, A., Estiar, M. A., Hashemkhani, M., Kiani, P., Mohamed, E.,
Saeidi, A., Taheri, M., Avci, P., et al. (2015). Carbon nanotubes part ii: a remarkable carrier
for drug and gene delivery. Expert opinion on drug delivery, 12(7):1089{1105.
Kaur, J., Gill, G. S., and Jeet, K. (2019). Applications of carbon nanotubes in drug delivery:
a comprehensive review. In Characterization and biology of nanomaterials for drug delivery,
pages 113{135. Elsevier.
Kauzmann, W. and Eisenberg, D. (1969). The structure and properties of water. Clarendon Press.
Kawasaki, T. and Kim, K. (2017). Identifying time scales for violation/preservation of stokes-
einstein relation in supercooled water. Science advances, 3(8):e1700399.
Kegel, W. K. and van Blaaderen, A. (2000). Direct observation of dynamical heterogeneities in
colloidal hard-sphere suspensions. Science, 287(5451):290{293.
97
Kell, G. (1972). Thermodynamic and transport properties of
uid water. In The Physics and
Physical Chemistry of Water, pages 363{412. Springer.
Kell, G. S. (1975). Density, thermal expansivity, and compressibility of liquid water from 0.
deg. to 150. deg.. correlations and tables for atmospheric pressure and saturation reviewed and
expressed on 1968 temperature scale. Journal of Chemical and Engineering Data, 20(1):97{105.
Kennedy, S. J. and Wheeler, J. C. (1983). On the density anomaly in sulfur at the polymerization
transition. The Journal of Chemical Physics, 78(3):1523{1527.
Khademi, M., Kalia, R. K., and Sahimi, M. (2015). Dynamics of supercooled water in nanotubes:
Cage correlation function and diusion coecient. Physical Review E, 92(3).
Khademi, M. and Sahimi, M. (2011). Molecular dynamics simulation of pressure-driven water

ow in silicon-carbide nanotubes. Journal of Chemical Physics, 135(20):7.
Khademi, M. and Sahimi, M. (2016). Static and dynamic properties of supercooled water in small
nanotubes. The Journal of Chemical Physics, 145(2):024502.
Khataee, A., Azamat, J., and Bayat, G. (2016). Separation of nitrate ion from water using silicon
carbide nanotubes as a membrane: insights from molecular dynamics simulation. Computational
Materials Science, 119:74{81.
Kheirabadi, A. M., Moosavi, A., and Akbarzadeh, A. M. (2014). Nano
uidic transport inside
carbon nanotubes. Journal of Physics D-Applied Physics, 47(6).
Kingston, C., Zepp, R., Andrady, A., Boverhof, D., Fehir, R., Hawkins, D., Roberts, J., Sayre,
P., Shelton, B., Sultan, Y., et al. (2014). Release characteristics of selected carbon nanotube
polymer composites. Carbon, 68:33{57.
Kirkwood, J. G. (1939). The dielectric polarization of polar liquids. The Journal of Chemical
Physics, 7(10):911{919.
Koga, K., Gao, G. T., Tanaka, H., and Zeng, X. C. (2001). Formation of ordered ice nanotubes
inside carbon nanotubes. Nature, 412(6849):802{805.
K ohler, M. H. and Bordin, J. R. (2018). Surface, density, and temperature eects on the wa-
ter diusion and structure inside narrow nanotubes. The Journal of Physical Chemistry C,
122(12):6684{6690.
K ohler, M. H., Bordin, J. R., da Silva, L. B., and Barbosa, M. C. (2017). Breakdown of the
stokes{einstein water transport through narrow hydrophobic nanotubes. Physical Chemistry
Chemical Physics, 19(20):12921{12927.
Kohler, M. H. and da Silva, L. B. (2016). Size eects and the role of density on the viscosity of
water conned in carbon nanotubes. Chemical Physics Letters, 645:38{41.
K ohler, M. H. and Gavazzoni, C. (2019). Water freezing in mos2 nanotubes. The Journal of
Physical Chemistry C, 123(22):13968{13975.
Kohlrausch, R. (1854). Theorie des elektrischen r uckstandes in der leidener
asche. Annalen der
Physik, 167(2):179{214.
Kolesnikov, A. I., Zanotti, J.-M., Loong, C.-K., Thiyagarajan, P., Moravsky, A. P., Loutfy, R. O.,
and Burnham, C. J. (2004). Anomalously soft dynamics of water in a nanotube: a revelation
of nanoscale connement. Physical review letters, 93(3):035503.
98
Kotzar, G., Freas, M., Abel, P., Fleischman, A., Roy, S., Zorman, C., Moran, J. M., and Melzak,
J. (2002). Evaluation of mems materials of construction for implantable medical devices. Bio-
materials, 23(13):2737{2750.
Kroto, H. W., Heath, J. R., O'Brien, S. C., Curl, R. F., and Smalley, R. E. (1985). C60:
Buckminsterfullerene. Nature, 318(6042):162.
Kumar, P., Buldyrev, S. V., Starr, F. W., Giovambattista, N., and Stanley, H. E. (2005). Ther-
modynamics, structure, and dynamics of water conned between hydrophobic plates. Physical
Review E, 72(5):051503.
Kyakuno, H., Matsuda, K., Yahiro, H., Inami, Y., Fukuoka, T., Miyata, Y., Yanagi, K., Maniwa,
Y., Kataura, H., Saito, T., Yumura, M., and Iijima, S. (2011). Conned water inside single-
walled carbon nanotubes: Global phase diagram and eect of nite length. Journal of Chemical
Physics, 134(24):14.
Laksmono, H., McQueen, T. A., Sellberg, J. A., Loh, N. D., Huang, C., Schlesinger, D., Sierra,
R. G., Hampton, C. Y., Nordlund, D., Beye, M., et al. (2015). Anomalous behavior of the
homogeneous ice nucleation rate in \no-man's land". The journal of physical chemistry letters,
6(14):2826{2832.
Larnaudie, S. C., Brendel, J. C., Romero-Canel on, I., Sanchez-Cano, C., Catrouillet, S., Sanchis,
J., Coverdale, J. P., Song, J.-I., Habtemariam, A., Sadler, P. J., et al. (2018). Cyclic peptide{
polymer nanotubes as ecient and highly potent drug delivery systems for organometallic
anticancer complexes. Biomacromolecules, 19(1):239{247.
Lau, K.-t., Gu, C., and Hui, D. (2006). A critical review on nanotube and nanotube/nanoclay
related polymer composite materials. Composites Part B: Engineering, 37(6):425{436.
Lee, W. J., Maiti, U. N., Lee, J. M., Lim, J., Han, T. H., and Kim, S. O. (2014). Nitrogen-doped
carbon nanotubes and graphene composite structures for energy and catalytic applications.
Chemical Communications, 50(52):6818{6830.
Liu, B., Wu, R., Baimova, J. A., Wu, H., Law, A. W.-K., Dmitriev, S. V., and Zhou, K. (2016a).
Molecular dynamics study of pressure-driven water transport through graphene bilayers. Phys-
ical Chemistry Chemical Physics, 18(3):1886{1896.
Liu, B., Wu, R., Law, A. W.-K., Feng, X.-Q., Bai, L., and Zhou, K. (2016b). Channel morphology
eect on water transport through graphene bilayers. Scientic reports, 6:38583.
Liu, L., Zhang, L., Sun, Z., and Xi, G. (2012). Graphene nanoribbon-guided
uid channel: a fast
transporter of nano
uids. Nanoscale, 4(20):6279{6283.
Liu, Y. and Wang, Q. (2005). Transport behavior of water conned in carbon nanotubes. Physical
review B, 72(8):085420.
Liu, Y., Wang, Q., Wu, T., and Zhang, L. (2005). Fluid structure and transport properties of
water inside carbon nanotubes. The Journal of chemical physics, 123(23):234701.
Lorenz, C. D., Chandross, M., Lane, J. M. D., and Grest, G. S. (2010). Nanotribology of water
conned between hydrophilic alkylsilane self-assembled monolayers. Modelling and Simulation
in Materials Science and Engineering, 18(3):034005.
Louden, P. B. and Gezelter, J. D. (2018). Why is ice slippery? simulations of shear viscosity of
the quasi-liquid layer on ice. The Journal of Physical Chemistry Letters, 9(13):3686{3691.
99
Lubchenko, V. (2009). Shear thinning in deeply supercooled melts. Proceedings of the National
Academy of Sciences, 106(28):11506{11510.
Ma, J., Michaelides, A., Alfe, D., Schimka, L., Kresse, G., and Wang, E. (2011). Adsorption and
diusion of water on graphene from rst principles. Physical Review B, 84(3):033402.
Ma, X., Cambr e, S., Wenseleers, W., Doorn, S. K., and Htoon, H. (2017). Quasiphase transition in
a single le of water molecules encapsulated in (6, 5) carbon nanotubes observed by temperature-
dependent photoluminescence spectroscopy. Physical review letters, 118(2):027402.
Mahbubul, I., Khan, M. M. A., Ibrahim, N. I., Ali, H. M., Al-Sulaiman, F. A., and Saidur, R.
(2018). Carbon nanotube nano
uid in enhancing the eciency of evacuated tube solar collector.
Renewable energy, 121:36{44.
Mahdizadeh, S. J. and Goharshadi, E. K. (2013). Natural gas storage on silicon, carbon, and
silicon carbide nanotubes: a combined quantum mechanics and grand canonical monte carlo
simulation study. Journal of Nanoparticle Research, 15(1):20.
Mahoney, M. W. and Jorgensen, W. L. (2000). A ve-site model for liquid water and the repro-
duction of the density anomaly by rigid, nonpolarizable potential functions. The Journal of
chemical physics, 112(20):8910{8922.
Malek, K. and Sahimi, M. (2010). Molecular dynamics simulations of adsorption and diusion of
gases in silicon-carbide nanotubes. Journal of Chemical Physics, 132(1).
Mallamace, F., Branca, C., Broccio, M., Corsaro, C., Mou, C. Y., and Chen, S. H. (2007a). The
anomalous behavior of the density of water in the range 30k < t < 373k. Proceedings of the
National Academy of Sciences of the United States of America, 104(47):18387{18391.
Mallamace, F., Branca, C., Corsaro, C., Leone, N., Spooren, J., Stanley, H. E., and Chen, S.-H.
(2010). Dynamical crossover and breakdown of the stokes- einstein relation in conned water
and in methanol-diluted bulk water. The Journal of Physical Chemistry B, 114(5):1870{1878.
Mallamace, F., Broccio, M., Corsaro, C., Faraone, A., Majolino, D., Venuti, V., Liu, L., Mou,
C. Y., and Chen, S. H. (2007b). Evidence of the existence of the low-density liquid phase in
supercooled, conned water. Proceedings of the National Academy of Sciences of the United
States of America, 104(2):424{428.
Mamontov, E., Burnham, C., Chen, S.-H., Moravsky, A., Loong, C.-K., De Souza, N., and
Kolesnikov, A. (2006). Dynamics of water conned in single-and double-wall carbon nanotubes.
The Journal of chemical physics, 124(19):194703.
Mao, Z. and Sinnott, S. B. (2000). A computational study of molecular diusion and dynamic

ow through carbon nanotubes. The Journal of Physical Chemistry B, 104(19):4618{4624.
Marti, J. and Gordillo, M. (2001). Temperature eects on the static and dynamic properties of
liquid water inside nanotubes. Physical Review E, 64(2):021504.
Martincic, M. and Tobias, G. (2015). Filled carbon nanotubes in biomedical imaging and drug
delivery. Expert opinion on drug delivery, 12(4):563{581.
Mashl, R. J., Joseph, S., Aluru, N. R., and Jakobsson, E. (2003). Anomalously immobilized water:
A new water phase induced by connement in nanotubes. Nano Letters, 3(5):589{592.
Mavrandonakis, A., Froudakis, G. E., Schnell, M., and M uhlh auser, M. (2003). From pure carbon
to silicon- carbon nanotubes: an ab-initio study. Nano letters, 3(11):1481{1484.
100
Mehrjouei, E., Akbarzadeh, H., Shamichali, A. N., Abbaspour, M., Salemi, S., and Abdi, P. (2017).
Delivery of cisplatin anti-cancer drug from carbon, boron nitride, and silicon carbide nanotubes
forced by ag-nanowire: A comprehensive molecular dynamics study. Molecular Pharmaceutics,
14(7):2273{2284.
Mei, F., Zhou, X., Kou, J., Wu, F., Wang, C., and Lu, H. (2015). A transition between
bistable ice when coupling electric eld and nanoconnement. The Journal of Chemical Physics,
142(13):134704.
Memarian, F., Fereidoon, A., Khodaei, S., Mashhadzadeh, A. H., and Ganji, M. D. (2017).
Molecular dynamic study of mechanical properties of single/double wall sicnts: Consideration
temperature, diameter and interlayer distance. Vacuum, 139:93{100.
Mendon ca, B. H., de Freitas, D. N., K ohler, M. H., Batista, R. J., Barbosa, M. C., and de Oliveira,
A. B. (2019). Diusion behaviour of water conned in deformed carbon nanotubes. Physica A:
Statistical Mechanics and its Applications, 517:491{498.
Menon, M., Richter, E., Mavrandonakis, A., Froudakis, G., and Andriotis, A. N. (2004). Structure
and stability of sic nanotubes. Physical Review B, 69(11):4.
Mishra, A. K. and Ramaprabhu, S. (2010). Magnetite decorated multiwalled carbon nanotube
based supercapacitor for arsenic removal and desalination of seawater. The Journal of Physical
Chemistry C, 114(6):2583{2590.
Molinero, V. and Moore, E. B. (2009). Water modeled as an intermediate element between carbon
and silicon. The Journal of Physical Chemistry B, 113(13):4008{4016.
Mondal, S. and Bagchi, B. (2019). Water in carbon nanotubes: Pronounced anisotropy in dielectric
dispersion and its microscopic origin. The journal of physical chemistry letters, 10(20):6287{
6292.
Mubarak, N., Sahu, J., Abdullah, E., and Jayakumar, N. (2014). Removal of heavy metals from
wastewater using carbon nanotubes. Separation & Purication Reviews, 43(4):311{338.
Mukherjee, B., Maiti, P. K., Dasgupta, C., and Sood, A. (2007). Strong correlations and ckian
water diusion in narrow carbon nanotubes. The Journal of chemical physics, 126(12):124704.
Nakaoka, S., Yamaguchi, Y., Omori, T., Kagawa, M., Nakajima, T., and Fujimura, H. (2015).
Molecular dynamics analysis of the velocity slip of a water and methanol liquid mixture. Physical
Review E, 92(2):022402.
Nanok, T., Artrith, N., Pantu, P., Bopp, P. A., and Limtrakul, J. (2009). Structure and dynamics
of water conned in single-wall nanotubes. The Journal of Physical Chemistry A, 113(10):2103{
2108.
Naserifar, S., Brooks, D. J., Goddard III, W. A., and Cvicek, V. (2017). Polarizable charge
equilibration model for predicting accurate electrostatic interactions in molecules and solids.
The Journal of chemical physics, 146(12):124117.
Naserifar, S. and Goddard, W. A. (2019). Liquid water is a dynamic polydisperse branched
polymer. Proceedings of the National Academy of Sciences, 116(6):1998{2003.
Naserifar, S. and Goddard III, W. A. (2018). The quantum mechanics-based polarizable force
eld for water simulations. The Journal of chemical physics, 149(17):174502.
101
Naserifar, S. and Goddard III, W. A. (2019). Anomalies in supercooled water at 230 k arise from
a 1d polymer to 2d network topological transformation. The Journal of Physical Chemistry
Letters, 10(20):6267{6273.
Naserifar, S., Oppenheim, J. J., Yang, H., Zhou, T., Zybin, S., Rizk, M., and Goddard III, W. A.
(2019). Accurate non-bonded potentials based on periodic quantum mechanics calculations
for use in molecular simulations of materials and systems. The Journal of chemical physics,
151(15):154111.
Nomura, K., Kaneko, T., Bai, J., Francisco, J. S., Yasuoka, K., and Zeng, X. C. (2017). Evidence
of low-density and high-density liquid phases and isochore end point for water conned to
carbon nanotube. Proceedings of the National Academy of Sciences, page 201701609.
Odom, T. W., Huang, J.-L., Kim, P., and Lieber, C. M. (1998). Atomic structure and electronic
properties of single-walled carbon nanotubes. Nature, 391(6662):62.
Ohba, T. (2014). Size-dependent water structures in carbon nanotubes. Angewandte Chemie
International Edition, 53(31):8032{8036.
Ou, B., Wang, J., Wu, Y., Zhao, S., and Wang, Z. (2019). Treatment of polyaniline wastewater by
coupling of photoelectro-fenton and heterogeneous photocatalysis with black tio2 nanotubes.
ACS Omega, 4(6):9664{9672.
Pajzderska, A., Gonzalez, M. A., Mielcarek, J., and Wasicki, J. (2014). Water behavior in mcm-41
as a function of pore lling and temperature studied by nmr and molecular dynamics simula-
tions. The Journal of Physical Chemistry C, 118(41):23701{23710.
P erez-Hern andez, N., Luong, T. Q., Febles, M., Marco, C., Limbach, H.-H., Havenith, M., P erez,
C., Roux, M. V., P erez, R., and Mart n, J. D. (2012). The mobility of water molecules through
hydrated pores. The Journal of Physical Chemistry C, 116(17):9616{9630.
Petrushenko, I. K. and Petrushenko, K. B. (2015). Mechanical properties of carbon, silicon carbide,
and boron nitride nanotubes: eect of ionization. Monatshefte f ur Chemie-Chemical Monthly,
146(10):1603{1608.
Pettersen, E. F., Goddard, T. D., Huang, C. C., Couch, G. S., Greenblatt, D. M., Meng, E. C.,
and Ferrin, T. E. (2004). Ucsf chimera - a visualization system for exploratory research and
analysis. Journal of Computational Chemistry, 25(13):1605{1612.
Pettersson, L. G. M., Henchman, R. H., and Nilsson, A. (2016). Water - the most anomalous
liquid.
Plimpton, S. (1995). Fast parallel algorithms for short-range molecular-dynamics. Journal of
Computational Physics, 117(1):1{19.
Pokropivnyi, V. and Silenko, P. (2006). Silicon carbide nanotubes and nanotubular bers: Synthe-
sis, stability, structure, and classication. Theoretical and Experimental Chemistry, 42(1):3{15.
Pugliese, P., Conde, M., Rovere, M., and Gallo, P. (2017). Freezing temperatures, ice nanotubes
structures, and proton ordering of tip4p/ice water inside single wall carbon nanotubes. The
Journal of Physical Chemistry B, 121(45):10371{10381.
Qin, L.-C., Zhao, X., Hirahara, K., Miyamoto, Y., Ando, Y., and Iijima, S. (2000). Materials
science: The smallest carbon nanotube. Nature, 408(6808):50.
102
Rabani, E., Gezelter, J. D., and Berne, B. (1999). Direct observation of stretched-exponential re-
laxation in low-temperature lennard-jones systems using the cage correlation function. Physical
review letters, 82(18):3649.
Rabani, E., Gezelter, J. D., and Berne, B. J. (1997). Calculating the hopping rate for self-
diusion on rough potential energy surfaces: Cage correlations. Journal of Chemical Physics,
107(17):6867{6876.
Raczyn ski, P., G orny, K., Dendzik, Z., Samios, J., and Gburski, Z. (2019). Modeling the impact
of silicon-carbide nanotube on the phospholipid bilayer membrane: Study of nanoindentation
and removal processes via molecular dynamics simulation. The Journal of Physical Chemistry
C, 123(30):18726{18733.
Rahbari, A., Poursaeidesfahani, A., Torres-Knoop, A., Dubbeldam, D., and Vlugt, T. J. (2018).
Chemical potentials of water, methanol, carbon dioxide and hydrogen sulphide at low tempera-
tures using continuous fractional component gibbs ensemble monte carlo. Molecular Simulation,
44(5):405{414.
Raicu, V. and Feldman, Y. (2015). Dielectric relaxation in biological systems: Physical principles,
methods, and applications. Oxford University Press, USA.
Raju, M., Van Duin, A., and Ihme, M. (2018). Phase transitions of ordered ice in graphene
nanocapillaries and carbon nanotubes. Scientic reports, 8(1):1{11.
Ramazani, F. and Ebrahimi, F. (2016a). Uncertainties in the capillary lling of heterogeneous
water nanochannels. The Journal of Physical Chemistry C, 120(23):12871{12878.
Ramazani, F. and Ebrahimi, F. (2016b). Water imbibition into nonpolar nanotubes with extended
topological defects. Chemical Physics, 476:23{28.
Ramin, L. and Jabbarzadeh, A. (2013). Eect of water on structural and frictional properties of
self assembled monolayers. Langmuir, 29(44):13367{13378.
Ramos-Alvarado, B., Kumar, S., and Peterson, G. (2016). Hydrodynamic slip in silicon nanochan-
nels. Physical Review E, 93(3):033117.
Rein ten Wolde, P., Ruiz-Montero, M. J., and Frenkel, D. (1996). Numerical calculation of the
rate of crystal nucleation in a lennard-jones system at moderate undercooling. The Journal of
chemical physics, 104(24):9932{9947.
Rosenbloom, A. J., Sipe, D. M., Shishkin, Y., Ke, Y., Devaty, R. P., and Choyke, W. J. (2004).
Nanoporous sic: A candidate semi-permeable material for biomedical applications. Biomedical
Microdevices, 6(4):261{267.
Ryckaert, J. P., Ciccotti, G., and Berendsen, H. J. C. (1977). Numerical-integration of cartesian
equations of motion of a system with constraints - molecular-dynamics of n-alkanes. Journal of
Computational Physics, 23(3):327{341.
Sahimi, M. and Ebrahimi, F. (2019). Ecient transport between disjoint nanochannels by a water
bridge. Physical review letters, 122(21):214506.
Sakai, V. G., Alba-Simionesco, C., and Chen, S. H. (2011). Dynamics of Soft Matter: Neutron
Applications. Springer Science & Business Media.
Saleh, T. A. (2016). Nanocomposite of carbon nanotubes/silica nanoparticles and their use for
adsorption of pb (ii): from surface properties to sorption mechanism. Desalination and Water
Treatment, 57(23):10730{10744.
103
Sauer, G. and Borst, L. (1967). Lambda transition in liquid sulfur. Science, 158(3808):1567{1569.
Scal, L., Fraux, G., Boutin, A., and Coudert, F.-X. (2018). Structure and dynamics of water
conned in imogolite nanotubes. Langmuir, 34(23):6748{6756.
Schulz, R., Lindner, B., Petridis, L., and Smith, J. C. (2009). Scaling of multimillion-atom
biological molecular dynamics simulation on a petascale supercomputer. Journal of Chemical
Theory and Computation, 5(10):2798{2808.
Sellberg, J. A., Huang, C., McQueen, T. A., Loh, N., Laksmono, H., Schlesinger, D., Sierra,
R., Nordlund, D., Hampton, C., Starodub, D., et al. (2014). Ultrafast x-ray probing of water
structure below the homogeneous ice nucleation temperature. Nature, 510(7505):381.
Setoodeh, A. R., Jahanshahi, M., and Attariani, H. (2009). Atomistic simulations of the buckling
behavior of perfect and defective silicon carbide nanotubes. Computational Materials Science,
47(2):388{397.
Shaat, M. and Zheng, Y. (2019). Fluidity and phase transitions of water in hydrophobic and
hydrophilic nanotubes. Scientic reports, 9(1):1{12.
Shank, A., Wang, Y., Kaledin, A., Braams, B. J., and Bowman, J. M. (2009). Accurate ab
initio and \hybrid" potential energy surfaces, intramolecular vibrational energies, and classical
ir spectrum of the water dimer. The Journal of chemical physics, 130(14):144314.
Shen, D., Zhan, Z., Liu, Z., Cao, Y., Zhou, L., Liu, Y., Dai, W., Nishimura, K., Li, C., Lin, C.-T.,
et al. (2017). Enhanced thermal conductivity of epoxy composites lled with silicon carbide
nanowires. Scientic Reports, 7(1):2606.
Shen, H. J. (2009). Thermal-conductivity and tensile-properties of bn, sic and ge nanotubes.
Computational Materials Science, 47(1):220{224.
Shi, R., Russo, J., and Tanaka, H. (2018). Origin of the emergent fragile-to-strong transition in
supercooled water. Proceedings of the National Academy of Sciences, 115(38):9444{9449.
Shimizu, H., Toyoda, K., Mawatari, K., Terabe, S., and Kitamori, T. (2019). Femtoliter gradient
elution system for liquid chromatography utilizing extended nano
uidics. Analytical chemistry,
91(4):3009{3014.
Singh, R. S., Biddle, J. W., Debenedetti, P. G., and Anisimov, M. A. (2016). Two-state ther-
modynamics and the possibility of a liquid-liquid phase transition in supercooled tip4p/2005
water. Journal of Chemical Physics, 144(14).
Soper, A. K. and Ricci, M. A. (2000). Structures of high-density and low-density water. Physical
Review Letters, 84(13):2881{2884.
Steinhardt, P. J., Nelson, D. R., and Ronchetti, M. (1983). Bond-orientational order in liquids
and glasses. Physical Review B, 28(2):784.
Striolo, A. (2006). The mechanism of water diusion in narrow carbon nanotubes. Nano letters,
6(4):633{639.
Swope, W. C., Andersen, H. C., Berens, P. H., and Wilson, K. R. (1982). A computer simulation
method for the calculation of equilibrium constants for the formation of physical clusters of
molecules: Application to small water clusters. The Journal of Chemical Physics, 76(1):637{
649.
104
Taghavi, F., Javadian, S., and Hashemianzadeh, S. M. (2013). Molecular dynamics simulation of
single-walled silicon carbide nanotubes immersed in water. Journal of Molecular Graphics and
Modelling, 44:33{43.
Takaiwa, D., Hatano, I., Koga, K., and Tanaka, H. (2008). Phase diagram of water in carbon
nanotubes. Proceedings of the National Academy of Sciences of the United States of America,
105(1):39{43.
Talla, J. A. (2019). Electronic properties of silicon carbide nanotube with stone wales defects
under uniaxial pressure: a computational study. Computational Condensed Matter, 19:e00378.
Teixeira, J. (2012). Recent experimental aspects of the structure and dynamics of liquid and
supercooled water. Molecular Physics, 110(5):249{258.
Terso, J. (1988). New empirical approach for the structure and energy of covalent systems.
Physical Review B, 37(12):6991.
Terso, J. (1989). Modeling solid-state chemistry: Interatomic potentials for multicomponent
systems. Physical Review B, 39(8):5566{5568.
Thomas, J. A. and McGaughey, A. J. (2009). Water
ow in carbon nanotubes: transition to
subcontinuum transport. Physical Review Letters, 102(18):184502.
Togaya, M. (1997). Pressure dependences of the melting temperature of graphite and the electrical
resistivity of liquid carbon. Physical review letters, 79(13):2474.
Tomo, Y., Askounis, A., Ikuta, T., Takata, Y., Seane, K., and Takahashi, K. (2018). Superstable
ultrathin water lm conned in a hydrophilized carbon nanotube. Nano letters, 18(3):1869{
1874.
Treacy, M. J., Ebbesen, T. W., and Gibson, J. M. (1996). Exceptionally high young's modulus
observed for individual carbon nanotubes. Nature, 381(6584):678.
Tu, Q., Yang, Q., Wang, H., and Li, S. (2016). Rotating carbon nanotube membrane lter for
water desalination. Scientic reports, 6:26183.
Tunuguntla, R. H., Henley, R. Y., Yao, Y.-C., Pham, T. A., Wanunu, M., and Noy, A. (2017).
Enhanced water permeability and tunable ion selectivity in subnanometer carbon nanotube
porins. Science, 357(6353):792{796.
Van Der Spoel, D. and van Maaren, P. J. (2006). The origin of layer structure artifacts in
simulations of liquid water. Journal of chemical theory and computation, 2(1):1{11.
Van Duin, A. C., Dasgupta, S., Lorant, F., and Goddard, W. A. (2001). Reax: a reactive force
eld for hydrocarbons. The Journal of Physical Chemistry A, 105(41):9396{9409.
Varanasi, S. R., Subramanian, Y., and Bhatia, S. K. (2018). High interfacial barriers at narrow
carbon nanotube{water interfaces. Langmuir, 34(27):8099{8111.
Verlet, L. (1967). Computer" experiments" on classical
uids. i. thermodynamical properties of
lennard-jones molecules. Physical review, 159(1):98.
Wang, X., Fu, F., Peng, K., Huang, Q., Li, W., Chen, X., and Yang, Z. (2020). Understanding of
competitive hydrogen bond behavior of imidazolium-based ionic liquid mixture around single-
walled carbon nanotubes. The Journal of Physical Chemistry C, 124(12):6634{6645.
105
Wang, X. and Liew, K. M. (2011). Hydrogen storage in silicon carbide nanotubes by lithium
doping. The Journal of Physical Chemistry C, 115(8):3491{3496.
Weeks, E. R., Crocker, J. C., Levitt, A. C., Schoeld, A., and Weitz, D. A. (2000). Three-
dimensional direct imaging of structural relaxation near the colloidal glass transition. Science,
287(5453):627{631.
Wei, N., Peng, X., and Xu, Z. (2014). Breakdown of fast water transport in graphene oxides.
Physical Review E, 89(1):012113.
Weiss, V. C., Rullich, M., K ohler, C., and Frauenheim, T. (2011). Kinetic aspects of the ther-
mostatted growth of ice from supercooled water in simulations. The Journal of chemical physics,
135(3):034701.
Wen, L., Li, F., and Cheng, H.-M. (2016). Carbon nanotubes and graphene for
exible electro-
chemical energy storage: from materials to devices. Advanced Materials, 28(22):4306{4337.
Whitby, M. and Quirke, N. (2007). Fluid
ow in carbon nanotubes and nanopipes. Nature
Nanotechnology, 2(2):87.
Widom, B. (1963). Some topics in the theory of
uids. The Journal of Chemical Physics,
39(11):2808{2812.
Williams, G. and Watts, D. C. (1970). Non-symmetrical dielectric relaxation behaviour arising
from a simple empirical decay function. Transactions of the Faraday society, 66:80{85.
Won, C. Y. and Aluru, N. (2008). Structure and dynamics of water conned in a boron nitride
nanotube. The Journal of Physical Chemistry C, 112(6):1812{1818.
Wu, R., Yang, M., Lu, Y., Feng, Y., Huang, Z., and Wu, Q. (2008). Silicon carbide nanotubes
as potential gas sensors for co and hcn detection. The Journal of Physical Chemistry C,
112(41):15985{15988.
Xu, B., Li, Y., Park, T., and Chen, X. (2011). Eect of wall roughness on
uid transport resistance
in nanopores. The Journal of chemical physics, 135(14):144703.
Xu, L., Mallamace, F., Yan, Z., Starr, F. W., Buldyrev, S. V., and Stanley, H. E. (2009). Ap-
pearance of a fractional stokes{einstein relation in water and a structural interpretation of its
onset. Nature Physics, 5(8):565{569.
Yang, L. and Guo, Y. (2020). Dynamics of water conned in a graphene nanochannel: dependence
of friction on graphene chirality. Nanotechnology, 31(23):235702.
Yang, R., Hilder, T. A., Chung, S.-H., and Rendell, A. (2011). First-principles study of water
conned in single-walled silicon carbide nanotubes. The Journal of Physical Chemistry C,
115(35):17255{17264.
Yeh, I.-C. and Hummer, G. (2004). System-size dependence of diusion coecients and viscosi-
ties from molecular dynamics simulations with periodic boundary conditions. The Journal of
Physical Chemistry B, 108(40):15873{15879.
Yonetani, Y. (2006). Liquid water simulation: A critical examination of cuto length. The Journal
of chemical physics, 124(20):204501.
Zang, J., Konduri, S., Nair, S., and Sholl, D. S. (2009). Self-diusion of water and simple alcohols
in single-walled aluminosilicate nanotubes. ACS nano, 3(6):1548{1556.
106
Zaragoza, A., Gonz alez, M. A., Joly, L., L opez-Montero, I., Canales, M., Benavides, A., and
Valeriani, C. (2019). Molecular dynamics study of nanoconned tip4p/2005 water: how con-
nement and temperature aect diusion and viscosity. Physical Chemistry Chemical Physics,
21(25):13653{13667.
Zeng, L., Zhou, X., Huang, X., and Lu, H. (2018). Phase transition-like behavior of the water
monolayer close to the polarized surface of a nanotube. Physical Chemistry Chemical Physics,
20(31):20391{20397.
Zhang, H. W., Ye, H. F., Zheng, Y. G., and Zhang, Z. Q. (2011). Prediction of the viscosity of
water conned in carbon nanotubes. Micro
uidics and Nano
uidics, 10(2):403{414.
Zhang, N., Song, Y., Huo, J., Li, Y., Liu, Z., Bao, J., Chen, S., Ruan, X., and He, G. (2017).
Formation mechanism of the spiral-like structure of a hydrogen bond network conned in a

uorinated nanochannel: A molecular dynamics simulation. The Journal of Physical Chemistry
C, 121(25):13840{13847.
Zhang, Y., Hu, K., Zhou, Y., Xia, Y., Yu, N., Wu, G., Zhu, Y., Wu, Y., and Huang, H. (2019). A
facile, one-step synthesis of silicon/silicon carbide/carbon nanotube nanocomposite as a cycling-
stable anode for lithium ion batteries. Nanomaterials, 9(11):1624.
Zhang, Y. F. and Huang, H. C. (2008). Stability of single-wall silicon carbide nanotubes - molecular
dynamics simulations. Computational Materials Science, 43(4):664{669.
Zhao, J. B., Qiao, Y., Culligan, P. J., and Chen, X. (2010). Conned liquid
ow in nanotube:
A numerical study and implications for energy absorption. Journal of Computational and
Theoretical Nanoscience, 7(2):379{387.
Zhao, K. and Wu, H. (2015). Fast water thermo-pumping
ow across nanotube membranes for
desalination. Nano Letters, 15(6):3664{3668.
Zhao, Z., Sun, C., and Zhou, R. (2020). Thermal conductivity of conned-water in graphene
nanochannels. International Journal of Heat and Mass Transfer, 152:119502.
Zheng, Y.-g., Ye, H.-f., Zhang, Z.-q., and Zhang, H.-w. (2012). Water diusion inside carbon
nanotubes: mutual eects of surface and connement. Physical Chemistry Chemical Physics,
14(2):964{971.
107

Abstract (if available)

Abstract Water under confinement is found everywhere: Inside fractured rock, within blood vessels and inside living tissue and cells. Nonetheless, the behavior of bulk and confined water still puzzles scientists around the world. There are many questions that have not been completely answered yet. ❧ The rise of the nanotechnology brought new materials with many potential applications. Carbon and Silicon Carbide Nanotubes (CNTs and SiCNTs) showed remarkable properties such as very high Young's modulus, extremely high thermal conductivity, capacity to be tuned and doped, and the potential to behave as metallic, non-metallic and semiconductors materials. Since the discovery of nanotubes two decades ago, researchers started the difficult task to study these materials and how the properties of water departs from the bulk when it is confined inside these nanotubes. ❧ Properties of water such as the melting and boiling points, viscosity, or even ice phases can suffer drastic changes when confined inside nanotubes. It is the purpose of this Thesis to identify and understand such changes. To achieve this goal we employed several computational techniques and models. Molecular Dynamics (MD) and Monte Carlo (MC) simulations were intensively used in this work. Similarly, accurate force fields such as the REBO potential and the RexPoN potential were also employed. ❧ Our methods and results are organized in this Thesis as follows: In Chapter 1, a brief explanation of the computational methods used in this Thesis is presented. Chapter 2 shows the results of our first work where we used water confined in an isolated CNT and showed that water can exhibit an ice-like behavior even when it reaches temperatures close to the bulk boiling point. In Chapter 3 we present how the confinement can change the rheology of water. We show that a non-Newtonian behavior appears in water when it is confined in a SiCNT. The non-Newtonian behavior appears ""early"" compared with bulk water. In Chapter 4 we seek to understand whether the changes in melting point when water is confined, are universal. We also studied the effect of the partial charges in the change of the melting point. Finally, in Chapter 5 we used the RexPoN potential to compute the excess chemical potential to use it in Monte Carlo simulations to create structures to be used in MD simulations. We found that there is a change in their melting point, but not as dramatic as initially though. The purpose of Chapter 5 is to construct a phase diagram using the RexPoN force field. ❧ We hope that the findings presented in this work will be useful to the scientific community and that nanofluidic systems can be developed based on our results.

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

Exploring properties of silicon-carbide nanotubes and their composites with polymers

PDF

Molecular dynamics studies of protein aggregation in unbounded and confined media

PDF

Hydrogen storage in carbon and silicon carbide nanotubes

PDF

Nanomaterials under extreme environments: a study of structural and dynamic properties using reactive molecular dynamics simulations

PDF

Light Emission from Carbon Nanotubes and Two-Dimensional Materials

PDF

Stability and folding rate of proteins and identification of their inhibitors

PDF

Large-scale molecular dynamics simulations of nano-structured materials

PDF

Machine-learning approaches for modeling of complex materials and media

PDF

A molecular dynamics study of interactions between the enzyme lysozyme and the photo-responsive surfactant azobenzene trimethylammonium bromide (azoTAB)

PDF

Molecular- and continuum-scale simulation of single- and two-phase flow of CO₂ and H₂O in mixed-layer clays and a heterogeneous sandstone

PDF

Chemical and mechanical deformation of porous media and materials during adsorption and fluid flow

PDF

Physics informed neural networks and electrostrictive cavitation in water

PDF

Printed and flexible carbon nanotube macroelectronics

PDF

Optoelectronic, thermoelectric, and photocatalytic properties of low dimensional materials

PDF

High throughput computational framework for synthesis and accelerated discovery of dielectric polymer materials using polarizable reactive molecular dynamics and graph neural networks

PDF

Atomistic modeling of the mechanical properties of metallic glasses, nanoglasses, and their composites

PDF

Quantum molecular dynamics and machine learning for non-equilibrium processes in quantum materials

PDF

Molecular-scale studies of mechanical phenomena at the interface between two solid surfaces: from high performance friction to superlubricity and flash heating

PDF

Molecular dynamics simulation study of initial protein unfolding induced by the photo-responsive surfactants, azoTAB

PDF

Efficient simulation of flow and transport in complex images of porous materials and media using curvelet transformation

Asset Metadata

Creator Cobeña-Reyes, José (author)

Core Title Dynamics of water in nanotubes: liquid below freezing point and ice-like near boiling point

School Viterbi School of Engineering

Degree Doctor of Philosophy

Degree Program Chemical Engineering

Publication Date 04/08/2021

Defense Date 03/18/2021

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag confined water,molecular simulation,nanotubes,OAI-PMH Harvest,thermodynamic properties

Language English

Contributor Electronically uploaded by the author (provenance)

Advisor Sahimi, Muhammad (committee chair), Kalia, Rajiv (committee member), Nakano, Aiichiro (committee member)

Creator Email jose.cobena.84@gmail.com

Permanent Link (DOI) https://doi.org/10.25549/usctheses-c89-439621

Unique identifier UC11666700

Identifier etd-CobenaReye-9427.pdf (filename),usctheses-c89-439621 (legacy record id)

Legacy Identifier etd-CobenaReye-9427.pdf

Dmrecord 439621

Document Type Dissertation

Rights Cobeña-Reyes, José

Type texts

Source University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Access Conditions The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the a...

Repository Name University of Southern California Digital Library

Repository Location USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA