Complex Behavior of Ordered and Icelike Water in Carbon Nanotubes

Aug 3, 2018 - The ordered water exists both below and above the bulk boiling point. The order is identified based on the ten Wolde parameter and may ...
0 downloads 0 Views 571KB Size
Subscriber access provided by Kaohsiung Medical University

Surfaces, Interfaces, and Catalysis; Physical Properties of Nanomaterials and Materials

Complex Behavior of Ordered and Icelike Water in Carbon Nanotubes Near the Bulk Boiling Point José Cobeña-Reyes, Rajiv K Kalia, and Muhammad Sahimi J. Phys. Chem. Lett., Just Accepted Manuscript • DOI: 10.1021/acs.jpclett.8b01953 • Publication Date (Web): 03 Aug 2018 Downloaded from http://pubs.acs.org on August 3, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 24 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry Letters

Complex Behavior of Ordered and Icelike Water in Carbon Nanotubes Near its Bulk Boiling Point Jos´ e Cobe˜ na-Reyes, Rajiv K. Kalia, and Muhammad Sahimi∗ Mork Family Department of Chemical Engineering and Materials Science, University of Southern California, Los Angeles, California, 90089-1211, USA

ABSTRACT: We report the results of extensive molecular dynamics (MD) simulation of water in a carbon nanotube (CNT) with a specific diameter over a wide range of temperature from 343 K to 423 K. In order to characterize the nature of water we have computed the Kirkwood g-factor, the ten Wolde parameter, the radial distribution, the cage correlation, and the intermediate scattering functions, the mean-square displacements of the water molecules, and the connectivity of the oxygen atoms. The computed properties provide evidence for a complex behavior. Some of the properties indicate an icelike structure, while others point to ordered (but not necessarily frozen) water. The connectivity of is close to 9. The ordered water exists both below and above its bulk boiling point. The order is identified based on the ten Wolde parameter and may explain, along with the dynamic slow down, the recent discovery of “ice” in CNTs near the bulk boiling point in a certain range of CNTs’ diameter, not seen in tubes of other sizes.

TOC Graphic

ACS Paragon Plus Environment

1

The Journal of Physical Chemistry Letters 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 24

It is well-known that some fluids exhibit anomalous behavior, which is different from both experimental observations and theoretical expectation for most liquids. For example, whereas the diffusion coefficient of most liquids decreases under compression, there are some fluids whose diffusivity increases under the same conditions. Examples include, but not limited to, Si,2,3 graphite,4 and liquid metals.5 The most important fluid with anomalous behavior is, however, water6 whose density is maximum at 277 K; it may possibly coexists within two distinct forms8 (see, however, Selberg et al.9 and Laksmono et al.10 ); it is in a supercooled state at sufficiently low temperatures,8 and the behavior of its isothermal compressibility, isobaric heat capacity and thermal expansion coefficient11 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 amorphous ice, and that, unlike most other liquids, volumetric expansion of water at 277 K and lower temperatures is due to the reduction in its entropy caused by tetrahedral symmetry of the local order around each water molecule12 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 confined media, and in particular in nanotubes? One cannot overstate the importance of the question, because understanding 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,11,12 as well as computer simulations indicate that the water flux in carbon nanotubes13 (CNTs) and silicon-carbide nanotubes14 (SiCNTs) may be three to four orders of magnitude larger than what is predicted by the classical Hagen-Poiseuille flow (HP) in cylindrical tubes. In addition, it has been demonstrated that other factors, such as functionalization of nanotubes, strongly influence water flow in CNTs.15,16 Similarly, whereas diffusion of water in unbounded media follows the Einstein relation and is affected by only the local temperature and pressure, in a confined medium the van der Waals and Coulombic interactions between water molecules and the medium’s walls influence their mobility, which typically reduces the rates of diffusion.17,18 The Stokes-Einstein relation between the diffusivity and viscosity of water also breaks down in confined media,19 including CNTs and SiCNTs.20,21 Using Raman spectroscopy, Agrawal et al.22 studied the behavior of water in isolated CNTs. They reported very high sensitivity to the diameter of CNTs and much larger increases in the ACS Paragon Plus Environment

2

Page 3 of 24 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry Letters

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.22 reported various unusual phenomena that depended on the size of CNTs, including reversible melting bracketed to 105151 ◦ C and 87-117 ◦ C, near-ambient phase transitions bracketed between 15-49 ◦ C and 3-30 ◦

C, and depression of the freezing point between −35 ◦ C and 10 ◦ C. Such an unusual behavior

had never been reported before. In this Letter 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.22 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.22 Results. We have computed several quantities i8n order to characterize the state of water in the CNTs. In what follows we present and discuss the results. The Kirkwood g-Factor. We computed the Kirkwood g-factor for two purposes. One was, as described in the Supplementary Information (SI), 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 for this purpose are presented in the SI. 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 1 presents the results. The g-factor at high temperatures T indicates no particular dipole orientation (see also Figure S2 in the SI). As T decreases, however, the g-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 temperatures, and the higher temperatures that begin at at 373 K, indicating different 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. ACS Paragon Plus Environment

3

The Journal of Physical Chemistry Letters 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 24

The ten Wolde Parameter. As explained in the section Methods, we computed the ten Wolde parameter dl (i, j) for l = 6. Figure 2 compares the results for water in CNTs at two temperatures with the corresponding results for bulk water. There are very significant differences between the confined water order parameter and that in the bulk. Note that values of d6 (i, j) 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. Dependence of the Potential Energy on the Radial Distance. An interesting feature is the dependence of the potential energy U (r) on the radial distance r. As Figure S3 in the SI indicates, 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 U (r) of each group. As Figure 3 indicates, the water molecules near the CNT’s wall possess much lower values of U (r) than those near the center, hence indicating a slow down in the dynamics in that layer. Such a difference in the dynamics was also confirmed by calculating the cage correlation (CC) function, which will be discussed next. The Cage Correlation Function. Figure 4(a) presents the computed CC functions C(t) at four temperatures. Although C(t) decays with time, but even after 40 ns it is still significantly 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 4(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 low potential energies than those near the ACS Paragon Plus Environment

4

Page 5 of 24 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry Letters

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 different phases in a CNT.1 Intermediate Scattering Function. Figure 5 presents Fs (κ, t) 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 5 are similar to those reported in several studies of supercooled or glasslike water24−26 and confined water at very low temperatures27,28 in various materials. The shape of Fs (κ, t) and the relaxation times that we have computed are also comparable to those in the aforementioned studies. As Figure 5 indicates, Fs (κ, t) behaves qualitatively similar to the CC function C(t) in that, it indicates a notable slowing down of the dynamics. That is, the water molecules go thorough two stages of relaxation. The first 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,25 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. (14) yielded, τl ≈ 88 ns for κ = 2, the highest peak in the structure factor. 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 6, we see that at T = 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 are around 10 − 12 ˚ A2 , 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 significant 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. Radial Distribution Function and Mean Connectivity of the Atoms. Figure 7 presents ACS Paragon Plus Environment

5

The Journal of Physical Chemistry Letters 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 24

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 differences. The main peak in g(r)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 about 3.6 ˚ A, which corresponds to the first well in g(r)OO . g(r)OH , on the other hand, contains information about the hydrogen bonds and, of course, its first peak corresponds to the O-H bond in a water molecule, so that the second peak indicates the most likely place to find a H atom from an O atom. This distance is about 1.73 ˚ A. On the other hand, looking at g(r)HH , we see that the first 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 connectivity Z 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. This paper 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.. We are currently carrying out extensive MD simulations in the CNTs with other sizes that were mentioned in the Introduction, in order to identify the type of ice that may form in them. We are also studying the same problems in two other types of nanotubes, namely, silicon-carbide and boron nitride nanotubes in order to address an important question: is icelike structure of water at such elevated temperatures in CNTs a universal phenomenon occurring in other types of nanotubes? Methods: Molecular Models and MD Simulation Protocol. The MD simulations were ACS Paragon Plus Environment

6

Page 7 of 24 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry Letters

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.22 The CNT was then inserted inside a simulation cell of size ˚3 , filled with enough water molecules to have a water density of 1 g/cm3 . We then 3502 × ` A 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 finally 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 (N V E) 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 simulations were carried out using the velocity-Verlet algorithm in a (N V T ) ensemble for 40 ns. All the simulations were carried out using the LAMMPS package,29 and the system was visualized using Chimera.30 The CNT contained 3620 carbon atoms, and there were 851 water molecules in the nanotube. The time step was always 1 fs. The REBO potential31 was used to model the CNT, while the TIP4P-EW32 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.32 The van der Waals interactions between the carbon atoms and the water molecules were represented by the LJ potential with a cutoff distance of 11 ˚ A. The LJ parameters for graphite33 , σ = 3.4 ˚ A and /kB = 28.1 K−1 , were utilized for the carbon atoms, where kB is the Boltzmann’s factor. √ The Lorentz-Berthelot mixing rules, σij = (σii + σjj )/2 and ij = ii jj , were utilized for the pairs ij. The Coulombic interactions were computed using the particle-particle-particle-mesh method. Computing the Kirkwood g-Factor. When simulating any phenomenon in strong confinement, 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 artificial crystallinity that may arise due to the periodicity and summations methods for ACS Paragon Plus Environment

7

The Journal of Physical Chemistry Letters 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 24

long-range interactions, such as the Ewald or particle-particle-particle-mesh method.34 Previous studies35−38 showed that the Kirkwood g-factor gK is an accurate tool of checking the existence of this type of artifact. gK is a measure of the orientational order of the dipole moment of the molecules, and represents the correlation of dipole µi of water molecule i with the total dipole of the water molecules located within a distance r from the molecule.37 gK = 1 when the dipole moment of the molecules lacks a specific orientation; for gK > 1 the dipole moments are parallel to an imposed field, whereas gK < 1 represents a case in which the dipole moments are aligned perpendicular to the field.39 In addition to gaining insight into the problem that we study, we also used gK to identify the right size of the simulation box. The g-factor is defined as39,40 gK (r) =

X µi · µj rij