Jump Reorientation of Water Molecules Confined in Narrow Carbon

Jul 8, 2009 - The angle of rotation of water molecule H1*O*H2* (see Figure 1a) was ...... Best , R. B.; Hummer , G. Reaction Coordinates and Rates fro...
15 downloads 0 Views 2MB Size
Download PDF
10322

J. Phys. Chem. B 2009, 113, 10322–10330

Jump Reorientation of Water Molecules Confined in Narrow Carbon Nanotubes Biswaroop Mukherjee,*,†,‡ Prabal K. Maiti,*,† Chandan Dasgupta,†,‡ and A. K. Sood§ Centre for Condensed Matter Theory, Department of Physics, Indian Institute of Science, Bangalore 560 012, India, Department of Physics, Indian Institute of Science, Bangalore 560 012, India, and Condensed Matter Theory Unit, Jawaharlal Nehru Centre for AdVanced Scientific Research, Bangalore 560 064, India ReceiVed: May 3, 2009

We used molecular dynamics (MD) simulations to study the reorientational dynamics of water molecules confined inside narrow carbon nanotubes immersed in a bath of water. Our simulations show that the confined water molecules exhibit bistability in their reorientational relaxation, which proceeds by angular jumps between the two stable states. The angular jump of a water molecule in the bulk involves the breaking of a hydrogen bond with one of its neighbors and the formation of a hydrogen bond with a different neighbor. In contrast, the angular jump of a confined water molecule corresponds to an interchange of the two hydrogen atoms that can form a hydrogen bond with the same neighbor. The free energy barrier between these two states is a few kBT. The analytic solution of a simplified two-state jump model that qualitatively explains the reorientational behavior observed in simulations is also presented. 1. Introduction Liquid water is a fascinating complex liquid, most of whose anomalous properties arise from its ability to form a hydrogenbonded three-dimensional network.1 This network of hydrogen bonds in liquid water is disordered and shows temporal variations primarily as a result of the rotation of the water molecules. The important processes of proton transfer and proton conduction2 through the hydrogen-bonded network and hydration of proteins and micelles3,4 depend very sensitively on the dynamics of these hydrogen bonds. It is believed that the ratelimiting step for the conduction of protons is the temperatureinduced breakage of hydrogen bonds between water molecules.2 Among the early efforts at understanding the hydrogen-bond structure of liquid water, the flickering cluster model of Frank et al.,5 the site percolation model of Stanley et al.,6 and the random network model of Rice et al.7 are notable. Although these models were able to explain some of the anomalies obseved in bulk water, it was left to molecular dynamics (MD) simulations by Rahman et al.8 to provide a more microscopic picture of some aspects of the statics and dynamics of liquid water. More recently, Luzar et al.9,10 tried to understand the temporal evolution of the hydrogen-bond population in terms of a diffusion model by matching the predictions of the model with those of classical MD simulations. Despite many years of numerical investigation, a clear picture of the molecular mechanism of the reorientation of water molecules has been lacking. Recently, Laage et al.11 showed, using classical MD simulations, that reorientational relaxation in bulk water at room temperature proceeds by large-amplitude angular steps, not by small diffusive steps. Laage et al.11-13 used a jump model developed by Ivanov14 to calculate the orientational relaxation times and compared these relaxation times to those obtained independently from MD simulations and experiments. The * To whom correspondence should be addressed. E-mail: biswa@ physics.iisc.ernet.in (B.M.), [email protected] (P.K.M.). † Centre for Condensed Matter Theory, Department of Physics, Indian Institute of Science. ‡ Jawaharlal Nehru Centre for Advanced Scientific Research. § Department of Physics, Indian Institute of Science.

Ivanov model14 requires two parameters, the amplitude of the jump and the mean waiting time between successive jumps. The values of these two parameters can be obtained from the MD simulations. Relaxation times obtained from the extended jump model (a version of the jump model that takes into account the relaxation of the oxygen atoms between jumps) agree well with those obtained from MD simulations. Progress in the field of infrared photon echo15 and time-resolved infrared spectroscopy now allows studies of the hydrogen-bond dynamics on the femtosecond time scale,16-19 which can serve as a testing ground of these ideas. It is known that confinement leads to significant changes in the structure and dynamics of water,20 owing to the fact that confinement modifies the hydrogen-bond structure found in bulk water. Recent experiments on the reorientational dynamics of water in reverse micelles [sodium bis(2-ethylhexyl) sulfosuccinate, or AOT] by spectrally resolved ultrafast spectroscopy21 yielded a wealth of information about the reorientational dynamics of confined water. Therefore, it is interesting to investigate how the mechanism of reorientational relaxation, specifically the jumps, is modified due to confinement. We studied the reorientational dynamics of water molecules inside very narrow carbon nanotubes. The narrow diameter of these nanotubes allows only a single file of water molecules inside. These water molecules form an orientationally ordered chain with the average dipole moment of the water molecules pointing either parallel or antiparallel to the axis of the nanotube.22,23 Each water molecule inside the chain is hydrogenbonded to two water molecules, one in front and one behind. This leaves one hydrogen atom for each molecule free, which does not participate in any hydrogen bond. On a longer time scale of several nanoseconds, the aligned dipoles relax by migration of hydrogen-bonding defects across the polarized water chain.24 These features make the reorientational dynamics of the confined water molecules significantly different from that in the bulk. The objective of our study was to explore the differences between the reorientational relaxation of the bulk and confined water molecules at shorter time scales. Our main results are summarized below.

10.1021/jp904099f CCC: $40.75  2009 American Chemical Society Published on Web 07/08/2009

Reorientation of Water Molecules Confined in Nanotubes The confined water molecules relax by a mechanism in which the hydrogen atom participating in the hydrogen bond and the free hydrogen atom repeatedly exchange their positions. Each of these exchanges results in an angular jump (see Figure 2 below), and it was found that these jumps typically occur through angles of about 90°. To complete a jump, a water molecule has to cross a free energy barrier of height ∼2.5kBT (this is the value for TIP3P25 model of water) at T ) 300 K, most of which is provided by electrostatic energy. Thus, the system exhibits “bistability” in its reorientational dynamics, which proceeds via thermally activated jumps between two locally stable states. This bistable behavior is very different from the mechanism of jumps in the bulk, where thermally activated jumps make the system explore many locally stable minima. The angular jump of a water molecule in the bulk involves the breaking of a hydrogen bond with one of its neighbors and the formation of a hydrogen bond with a different neighbor. In contrast, the angular jump of a confined water molecule corresponds to an interchange of the two hydrogen atoms that each can form a hydrogen bond with the same neighbor. We also computed the distributions of the waiting time between two successive jumps and the time required to complete a successful jump. The distribution of the waiting time between successive jumps of the confined water molecules peaks at approximately 60 fs and exhibits a long, slowly decaying tail, with a mean waiting time of ∼1 ps. Our calculations of the reorientational correlation functions reveal features arising from the jump mechanism. To explain the peculiarities of the observed correlation functions, we also present a simplified “two-state” Kubo rotor model of this jump mechanism and calculate its correlation functions. The behavior of the correlation function is qualitatively similar to the results obtained from MD simulations. To the best of our knowledge, this is the first study of the details of the relaxation process via jumps for nanoconfined single-file water molecules. We performed these calculations for TIP3P,25 SPC/E,26 and TIP5P27 water models. This article discusses all of the results for the TIP3P water model in detail, and the final section summarizes the similarities and the differences between the results obtained using the other water models. The rest of this article is organized as follows: The details of the simulation and the methods used for analysis are described in section 2, the results obtained using the TIP3P water model are presented in section 3, section 4 discusses the similarities and differences between the results obtained using other water models, and section 5 contains our conclusions. 2. Simulation Details Open-ended (6,6) carbon nanotubes of various lengths were immersed in a bath of TIP3P25 water molecules. The simulations were performed using AMBER 7.28 During the simulations, the nanotubes were held fixed in a stiff harmonic external potential inside the simulation box, with their long axes along the z direction. Simulations were performed for nanotubes with a size of 24 unit cells (576 atoms, 56 Å). The simulation box contained 2223 water molecules. Initially, an NPT simulation was performed at a temperature of 300 K and 1 atm pressure for 10 ns with temperature regulation achieved using the Berendsen weak coupling method (0.5-ps time constant for heat bath coupling and 0.2-ps pressure relaxation time). Finally, NVE simulations were performed with a time step of 0.5 fs for a total time of 215 ps. Trajectories of the NVE simulations were used to calculate the geometrical and dynamical features of the jump. The interactions between various atoms were described

J. Phys. Chem. B, Vol. 113, No. 30, 2009 10323

Figure 1. (a) Definition of hydrogen bond between water molecules O*H1*H2* and O′H1′H2′. A hydrogen bond between these two water molecules is assumed to exist if the distance RO*O′ is less than 3.5 Å and if one of the angles H1*O*O′ or H2*O*O′ is less than 30°. (b) Definition of a correct jump: Between times t0 and t1, fH1(t) goes from 1 to 0 and back to 1, with fH2(t) remaining 0 throughout. This is not a jump. At time t2, fH1(t) goes from 1 to 0, and at a later instant t3, fH2(t) goes from 0 to 1. This is an example of a successful jump, starting at time t ) t2 and ending at time t ) t3.

by classical force fields; the carbon and oxygen atoms were modeled as Lennard-Jones particles with parameter values given in our earlier report.29 The carbon-carbon bond length was 1.4 Å, and the corresponding spring constant was 938 kcal/(mol Å2); the equilibrium C-C-C angle was 2π/3 rad, and the corresponding spring constant was 126 kcal/(mol rad2). The electrostatic interactions were calculated with the particle-mesh Ewald (PME) method using a cubic B-spline interpolation of order 4 and a 10-4 tolerance set for the direct-space sum cutoff. A real-space cutoff of 9 Å was used for both the electrostatic and van der Waals (vdW) interactions, with a nonbonded-list update frequency of 10. The initially empty nanotube was spontaneously filled within a few hundred picoseconds and remained filled with, on average, about 21 water molecules. The chain of water molecules was tightly packed inside the nanotube with an average density about 4 times that of bulk water.30 To characterize the jumps, we carried out the following procedure: A water molecule confined inside the nanotube was defined as being in either the “up” state or the “down” state, depending on the z component of the dipole moment of the particular water molecule. Suppose that the water molecules were in the up state. We searched the water molecules in ascending order of their axial coordinates, and for each water molecule at time instant t, we determined which hydrogen atom (H1 or H2) of the water molecule was participating in hydrogen bonding with the water molecule in front of it. If atom H1 was participating in hydrogen bonding, we put a flag, fH1(t), “on” [fH1(t) ) 1]; otherwise, fH1(t) ) 0. fH2(t) was defined in a similar way. We checked two geometrical conditions between two neighboring water molecules, H1*O*H2* and H1′O′H2′ (see Figure 1a). A hydrogen bond was assumed to exist if RO*O′ < 3.5 Å and θH*O*O′ < 30° (similar to the approach used in ref 11), where RO*O′ is the distance between the donor and the acceptor oxygen atoms and θH*O*O′ is the angle between the OH bond and the O*O′ vector. In this way, we obtained the functions fH1(t) and fH2(t) for each water molecule in the chain. A successful jump was defined to occur for a given molecule in the chain if its fH1(t) value went from 1 to 0 at some instant and its fH2(t) value went from 0 to 1 at a later instant (or vice versa). This process is schematically shown in Figure 1. Between times t0 and t1, fH1(t) goes from 1 to 0 and back to 1, with fH2(t) remaining 0 throughout. This is not a jump. At time t2, fH1(t) goes from 1 to 0, and at a later instant t3, fH2(t) goes from 0 to 1. This is an example of what we consider a successful jump. In this way, we obtained the starting time and the ending time of each successful jump, and from these data, the

10324

J. Phys. Chem. B, Vol. 113, No. 30, 2009

Figure 2. Large-amplitude angular jump of a confined water molecule. The molecule at the center (O*H1*H2*) is the one performing the large angular jump. It breaks one hydrogen bond at the beginning of the jump (A), when H2* is bound to oxygen O′, and forms another at the end (C), when H1* is bound to O′.

distributions of the waiting time between two jumps and the time required to complete a jump were computed. The angle of rotation of water molecule H1*O*H2* (see Figure 1a) was calculated in the following way: The angle that one of the OH vectors of the water molecule H1*O*H2* made with the line O*O′ at time instant t ) t2, which corresponds to the start of the jump, was calculated (θ2). The angle that the same OH vector made with the O*O′ line was also calculated at a later instant t ) t3 (θ3). The difference between these two angles, θ2 and θ3, was taken as the jump angle.

Mukherjee et al.

Figure 3. Temporal evolution of various distances and angles during a molecular jump. The plots in this figure were obtained from coherent averaging (performed by matching the time instant corresponding to the center of each jump) over about 100 successful jump trajectories. The trajectories chosen for the calculation of the average were those that did not have another jump within 250 fs before or after the time of the selected jump.

3. Results and Discussion Figure 2 shows three stages of a large-amplitude angular jump of the central, rotating water molecule O*H1*H2*. The hydrogen atom participating in the hydrogen bond changes identity in this process. The bistable nature of the orientational relaxation becomes clear from Figure 2. If the result of one jump is such that atom H1* is hydrogen-bonded and H2* is free, then the next jump must result in H2* being hydrogenbonded and H1* being free. The water molecule has to cross a free energy barrier of height ∼2.5kBT when it goes from one minimum to the other. The details of the free energy and the energy calculation are discussed below. In a total simulation time of 215 ps, with a time step of 0.5 fs, we detect about 2000 successful jumps, which were used in the calculation. The total time of 215 ps was divided into 5 parts, and for each portion, 16 water molecules that stayed inside the nanotube throughout this interval were tracked. Calculations were done on the trajectories of these 5 sets of 16 water molecules. Figure 3 shows the temporal evolution of the distances rH1*O′, rH2*O′, and rO*O′ (top panel) and the angle θO′O*H* (bottom panel) that one OH bond of the rotating water molecule makes with the line joining the atoms O* and O′ (see Figure 1a for the definitions) across the angular jump. From the top panel of Figure 3, it is clear that the distance between the oxygen atoms O* and O′ (rO*O′) remains almost constant, whereas the distances of H1* and H2* from O′ (rH1*O′ and rH2*O′) change across the jump. The bottom panel of Figure 3 shows how the angle θO′O*H* changes with time during the angular jump. A coherent averaging of all of the trajectories was done by matching the times corresponding to the centers of all of the jumps. For this plot, only those jumps that did not have any other jump within two time windows of 250 fs before and after this particular jump were chosen. About 100 such jumps were used to generate this average trajectory. The top panel in Figure 3 shows oscillations in the distance between one of the hydrogen atoms of the

Figure 4. One-dimensional free energy profile, ∆F ≈ -kBT ln[p(θO′O*H*)], calculated from the probability distribution of the angles O′O*H1* and O′O*H2*. The two stable states are near θO′O*H* ) 10° and θO′O*H* ) 90°, separated by a barrier at around 50°. The value of the barrier is about 2.5 kBT.

rotating water molecule (the one that does not form a hydrogen bond) and the oxygen atom of the water molecule in front. This is understandable: as this hydrogen atom is not bound, its fluctuations are expected to be greater than those of the bound hydrogen atom. The free energy barrier between these two locally stable states was estimated. The calculation was performed by evaluating the probability distribution of the angles O′O*H1* and O′O*H2* and then calculating the quantity

∆F ≈ -kBT ln[p(θO′O*H*)]

(1)

where T ) 300 K and p(θO′O*H*) is the probability distribution. The free energy profile is shown in Figure 4, which shows that the two stable states are near θO′O*H* ) 10° and θO′O*H* ) 90°, separated by a barrier at around 50°. The height of the barrier is about 2.5kBT. Comparing this result with the calculation of the energy barrier, which will be discussed below, we found that most of the contribution arises from the electrostatic energy, with the entropy contributing a small fraction. The value of the

Reorientation of Water Molecules Confined in Nanotubes

J. Phys. Chem. B, Vol. 113, No. 30, 2009 10325

Figure 5. Definition of the three relevant axes of rotation. Axes 1 and 3 are shown; axis 2 is perpendicular to the plane of the page. Figure 6. Energy as a function of the angle of rotation about axis 3.

free energy barrier in this case was less than the barrier for angular jumps calculated for bulk SPC/E water.31 To estimate the fraction of the barrier height contributed by the electrostatic energy, we performed a separate NVT simulation at 300 K for a system of 40 water molecules confined inside a (6,6) carbon nanotube of length 123 Å. Periodic boundary conditions were applied in the direction of the axis of the tube. The surrounding bulk water was absent in this case. We took a typical configuration from the simulation and minimized the energy. The resulting final structure was one in which most of the water molecules participated in two hydrogen bonds, donating one hydrogen bond to the water molecule in front of it and accepting one hydrogen bond from the water molecule behind it. We selected such a double-hydrogen-bonded water molecule and applied a set of rotations to it, keeping all of the other water molecules fixed. The axes of the various rotations always passed through the oxygen atom, ensuring that this atom did not undergo any displacement. This method ensured that vdW energy changes (which result from changes of the positions of the oxygen atoms only, as the hydrogen atoms in TIP3P water do not have any vdW contribution) were avoided and only the electrostatic component resulting from changes of the positions of the hydrogen atoms due to rotations was measured. Figure 5 defines the various axes of rotation used in our calculations. The axis marked 1 in Figure 5 is in the direction O*O′. Rotations about this axis keep hydrogen atom H1 hydrogen-bonded and atom H2 unbonded. Axis 2 is perpendicular to the plane of the page, and rotations about this axis result in transitions from a state where atom H1 is hydrogenbonded to a state where atom H2 is bonded and vice versa. Rotations about axis 3 take the plane of the water molecule away from the plane of the page and disrupt the hydrogen-bond configuration. Hence, a rotation about axis 3 is energetically expensive. To estimate this energy cost, we calculated the total potential energy of the rotating water molecule as a function of the angle of rotation about axis 3 as the angle of rotation varied from -π/2 to π/2. It is apparent from Figure 6 that rotations about axis 3 are energetically quite costly as they require the energy equivalent of about one hydrogen bond (∼7-8 kcal/ mol). To scan the energy as a function of the angles of rotation about axes 1 and 2, we used the following protocol: A rotation of the particular water molecule chosen for this calculation by an angle of π/2 about axis 2 resulted in a transition from a state where atom H1 was bonded (H2 free) to a state where atom H2 was bonded (H1 free). For each small rotation (the total rotation of angle π/2 about axis 2 was divided into 10 steps) about axis 2, we scanned all possible rotation angles (from 0 to 2π) about axis 1. Figure 7 shows the results of this exercise.

Figure 7. Two-dimensional potential energy surface of a rotating water molecule in the chain as a function of the two rotations about axes 1 and 2, which are defined in the text. The axis ∆θ1 goes from 0 to π, and ∆θ2 goes from 0 to π/2.

The starting point where ∆θ1 ) ∆θ2 ) 0 in Figure 7 represents the initial position of the water molecule that was chosen for the application of various rotations. It is evident from this figure that a rotation about axis 2 alone does not result in a transition from one energy minimum to another. This is clear if one observes the ∆θ1 ) 0 cut of the two-dimensional plot. One goes from the original energy minimum to higher-energy states as ∆θ2 is increased from 0 to π/2. If we also include rotations about axis 1 for each elementary angular displacement about axis 2, we see that there exists a low-energy path over a saddle point that connects the H1-bonded and H2-bonded states (Figure 2). Physically, this means that the plane of the water molecule also changes as it makes the jump. The effective energy barrier was calculated by following the energy that was the minimum of all the ∆θ1 configurations for each value of ∆θ2. The resulting barrier height is about 1.1 kcal/mol (∼2kBT at T ) 300 K). However, it should be pointed out that the details of the energy surface are sensitive to the position and orientation of the other water molecules in the water chain, which were kept fixed in the present method of calculation. The two locally stable states (see Figure 2) of each water molecule in the chain and the large-amplitude angular jumps between these two states strongly affect the reorientational dynamics of the confined water molecules. A qualitatively similar energy landscape is also exhibited by small water clusters

10326

J. Phys. Chem. B, Vol. 113, No. 30, 2009

Mukherjee et al.

〈τw〉 )

Figure 8. Cumulative distributions, w′(t) and j′(t), which have been defined in the text. j′(t) approaches unity by about 200 fs, whereas w′(t) exhibits a significant effect of the slowly decaying tail of the waiting-time distribution.

in an uniform electric field.32 The water molecule relaxes by relatively fast and small angular displacements in one of these two minima and also by jumps from one of these minima to the other, assisted by thermal fluctuations. This behavior is similar to that found in dense liquids where reorientational relaxation involves thermally activated transitions between different local minima of the energy or free-energy landscape, although the bistability discussed above is a feature special to this particular situation. Even in bulk liquids comprising orientable molecules, the local potential energy suface comprises local minima, but the positions of these minima change relatively quickly. In the present situation, the strong confinement forces the positions of the potential energy minima, which are primarily determined by the positions of the heavy oxygen atoms of the water molecules, to relax very slowly. The distribution of the waiting time between two successful jumps and the distribution of the time required to complete a jump were also computed. The waiting-time distribution peaks at about 60 fs and exhibits a long, slowly decaying tail that has a significant weight. The distribution of the time required to complete a jump peaks at 40 fs and does not have the long tail seen in the waiting-time distribution. To clearly demonstrate the importance of the tail of the waiting-time distribution we calculated and plotted the following quantities

w'(t) ) j'(t) )

∫0t w(t1) dt1 ∫0t j(t1) dt1

(2)

where w(t) and j(t) are the waiting-time and jump-time distributions, respectively. The quantity w′(t) is the area under the waiting-time distribution between times 0 and t; j′(t) is defined similarly. It is expected that w′(t) and j′(t) will start from 0 at small times and approach unity (for normalized distributions) at long times. Figure 8 shows w′(t) (dashed line) and j′(t) (solid line), and it is clear that j′(t) approaches unity much earlier (∼200 fs) than w′(t), which is very slow in approaching unity. This means that, even though the waiting-time distribution has a short-time peak at 60 fs, the effect of the slowly decaying tail of the distribution is significant. The mean waiting time, estimated from the formula

∫0∞ tw(t) dt ∫0∞ w(t) dt

(3)

is about 1 ps. This is roughly one-half of the mean waiting time between jumps (∼2 ps) in bulk SPC/E water. The decrease in the value of the mean waiting time is consistent with the fact that the free energy barrier for jumps under confinement is less than that for jumps in the bulk. We also calculated more conventional dynamical quantities such as the reorientational correlation function of the confined water molecules. The most widely accepted framework for understanding rotational relaxation originates from the work of Debye.33 The physical picture underlying this framework views rotational diffusion as a succession of small random steps. A unit vector fixed to the center of mass of a molecule is assumed to undergo random motion on an unit sphere. Solving the diffusion equation for the temporal evolution of the probability, P(φ,t), that a molecule reorients by an angle φ during time t yields the expression ∞

P(φ, t) )

+1 ∑ ( 2l4π )Pl[cos φ(t)] exp[-l(l + 1)Drt] l)0

(4) where Pl is the lth Legendre polynomial and Dr is the rotational diffusion constant in units of inverse time. When φ(t) is replaced by cos-1[eRi (t) · eRi (0)], the angular reorientation of a vector eRi fixed to the ith molecule, the lth-order correlation function is given by

〈∑ 〈∑ N

ClR(t) )

i)1 N

i)1

〉 〉

Pl[eiR(t) · eiR(0)]

(5)

Pl[eiR(0) · eiR(0)]

Of particular experimental interest is the correlation function corresponding to l ) 2, which can be probed by NMR experiments.34 In our calculations, the vector eR can be any vector joining two points on the water molecule. It is evident from eq 5 that the correlation time for the lth-order correlation function is τl ) [l(l + 1)Dr]-1 and the ratio of the relaxation times corresponding to l ) 1 (τ1) and l ) 2 (τ2) should be equal to 3, if the relaxation is Debye-like. Deviations from this relationship will occur if the relaxation involves nondiffusive processes. In particular, if large-amplitude angular jumps are present, then this ratio decreases from its ideal value of 3, as observed by Hynes et al.11 and Lombardo et al.35 The reason behind the deviation of this ratio from its value for diffusive orientation can be understood in the following way: Cl(t) decays when molecules have rotated sufficiently so that the appropriate lth-order spherical harmonic passes through its first zero. If the jumps are short, then it takes many jumps to reach the first zero, more for low-order spherical harmonics than for high-order ones, whereas if the jumps are long, then a small number of jumps carry the molecules through the first zero regardless of the order of the spherical harmonic. From Figure 9, it is evident that C1OH(t) does not decay to zero at long times, whereas C2OH(t) does. To understand the asymptotic values of C1OH(t) and C2OH(t), we present a small

Reorientation of Water Molecules Confined in Nanotubes

J. Phys. Chem. B, Vol. 113, No. 30, 2009 10327 constant value. To estimate this constant value, we observe that Figure 10 primarily comprises two peaks roughly at angles of 27° and 76°. If we assume that P(θ) ) 1/2[δ(θ - 27°) + δ(θ 76°)], then 〈cos θ〉2 ≈ 0.3, which is close to the saturation value of C1OH(t) in Figure 9. Similarly, the asymptotic value of the second-order correlation function can be determined from

COH 2 (t f ∞) )

OH Figure 9. Time evolution of COH 1 (t) and C2 (t). The fits of the correlation functions are between 0.6 and 4 ps. C1OH(t) ≈ 0.26 + -t/2 e-t/1.75, yielding τ1 as 1.75 ps. COH , yielding τ2 as 2 ps. The 2 (t) ≈ e two-state Kubo rotor model, presented in the text, qualitatively explains the features of the correlation functions obtained from MD simulations.

Figure 10. Distribution of the angles between the OH vector and the axis of the nanotube.

calculation. Let us denote the orientation of one OH vector of a water molecule at times 0 and t by n(0) ) [cos θ(0)]k + [sin θ(0) sin φ(0)]j + [sin θ(0) cos φ(0)]i and n(t) ) [cos θ(t)]k + [sin θ(t) sin φ(t)]j + [sin θ(t) cos φ(t)]i, respectively, where i, j, and k denote the laboratory-frame unit vectors. Hence, C1OH(t) is given by

COH 1 (t) ) 〈n(0) · n(t)〉 ) 〈cos θ(0) cos θ(t)〉 + 〈sin θ(0) sin θ(t)[cos φ(0) cos φ(t) + sin φ(0) sin φ(t)]〉 (6) At long times, the correlations cease; hence, 〈cos θ(0) cos θ(t)〉 can be written as 〈cos θ〉2, where the average is taken over stationary distribution of θ, where θ denotes the angle the OH vector makes with k. Figure 10 shows the time-averaged distribution of θ. The time-averaged distribution of φ for a particular value of the variable θ is uniform. This implies that 〈cos φ〉 and 〈sin φ〉 are 0; hence, the second term in the above equation drops out. This implies that the asymptotic value of the C1OH(t) is 〈cos θ〉2; that is, C1OH(t) decays not to 0, but to a

〈

〉

3[n(0) · n(t)]2 - 1 ) 2 3 3〈cos2 θ〉2 + 〈sin2 θ〉2 - 1 2 2

(7)

Substituting P(θ) ) 1/2[δ(θ - 27°) + δ(θ - 76°)], we obtain OH COH 2 (tf∞) ≈ 0.007. These observations suggest that C1 (t) does (t) does. To extract not decay to zero at long times, whereas COH 2 (t) was fitted to a constant and long-time-decay time scales, COH 1 an exponential, whereas C2OH(t) was fitted to an exponential. Figure 9 shows the result of the fits of the correlation functions between 0.6 and 4 ps. C1OH(t) ≈ 0.26 + e-t/1.75, yielding τ1 as 1.75 ps. The asymptotic value of 0.26 is quite close to our earlier estimate of 0.3. C2OH(t) ≈ e-t/2, yielding τ2 as 2 ps. The ratio of τ1 and τ2 in this case is about 0.88, which is much less than the diffusive-limit value of 3. The present observations are rationalized below, by deriving the correlation functions for a bistable Kubo rotor. One must note that this simplified model does not exactly mimic the situation of the MD simulations. Most importantly, it is a calculation involving a two-dimensional rotor, whereas the rotations of the water molecules occur in three dimensions. However, this model is presented to demonstrate the important effect of jumps, whose angles alternate in sign, which is an observation from the simulations. In this simplified model, we obtain correlations functions C1(t) and C2(t) and show that C1(t) does not decay to zero at long times, whereas C2(t) does, a behavior that is qualitatively similar to what we observed in our simulations. We consider the reorientational correlation functions for a two-state Kubo rotor,36 where the two states correspond to the situations when one OH bond is hydrogen-bonded to the water molecule in front and when it is free. From the calculation of the geometrical quantities averaged over trajectories, we know that the magnitude of the angular jump is approximately π/2. For a two-dimensional Kubo rotor, where each jump changes the angle by an amount ∆k, the correlation functions can be cast in the form

[∑ ∞

n

Cl(t) ) R

exp(il

n)0

]

∑ ∆k)P(n, t)

k)1

(8)

where P(n,t), is the probability that there are n jumps within a time interval t. P(n,t) can be assumed to have a Poisson distribution

P(n, t) )

(t/τ)n exp(-t/τ) n!

(9)

where τ is the average waiting time between two jumps. In the present problem of water molecules confined inside narrow nanotubes, the rotating vector corresponds to one of the OH bonds of a water molecule in the chain. The crucial feature

10328

J. Phys. Chem. B, Vol. 113, No. 30, 2009

Mukherjee et al.

special to the present case is the fact that the jumps of a particular water molecule are not independent but are strongly correlated. The angular jumps are, in fact, alternating in nature. If the first jump is through an angle of ∆, the next jump must be through an angle of -∆, the next one again through an angle of ∆, and so on. The equation for the correlation function will thus have the following special form in the present case

Cl(t) ) R[e0P(0, t) + eil∆P(1, t) + e0P(2, t) + ...]

(10) The terms involving ∆ and those without it can be collected separately, and the above equation takes the form

Cl(t) ) R{[P(0, t) + P(2, t) + ...] + eil∆[P(1, t) + P(3, t) + ...]} (11)

Figure 11. Temporal evolution of various distances and angles for TIP5P water inside a carbon nanotube.

It can be easily shown that

1 [P(0, t) + P(2, t) + ...] ) e-(t/τ) (e(t/τ) + e-(t/τ)) 2 1 [P(1, t) + P(3, t) + ...] ) e-(t/τ) (e(t/τ) - e-(t/τ)) (12) 2 Collecting all the terms, we can rewrite the correlation function as

[

Cl(t) ) R

(1 + e-2t/τ) (1 - e-2t/τ) + eil∆ 2 2

]

(13)

From previous calculations, we know that the magnitude of the jump is approximately π/2. Because cos(∆) ) 0 and cos(2∆) ) -1, we can write

1 e-2t/τ + 2 2 -2t/τ C2(t) ) e

C1(t) )

Figure 12. Temporal evolution of various distances and angles for SPC/E water inside a carbon nanotube.

TABLE 1: Comparison of Jump Parameters for Various Modelsa

(14)

This simple model calculation thus predicts that, at long times, C1(t) does not decay to 0, whereas C2(t) decays to 0 with a time scale of τ/2. The result that C1(t) does not decay to 0 is due to the fact that C1(t) is made up of a series of 1’s and 0’s: 1 when the initial and final positions of the vector coincide and 0 when these two positions are perpendicular. Hence, the average is 1 /2. Similarly, C2(t) goes to 0 at long times because it is made up of +1’s and -1’s. Other details of this reorientational behavior and the confinement-induced magnification of the relaxation anisotropy, which is not present in bulk water, are discussed in an earlier work.37 Experimentally, such reorientational anisotropies are also observed in the orientational dynamics of fluids such as carbon disulfide38 and methyl iodide39 confined inside nanoporous glasses, which have been studied using optical Kerr effect (OKE) spectroscopy. 4. Jumps in Confined TIP5P and SPC/E Water Whereas all of the results presented in this article thus far have been for the TIP3P model of water, similar analyses were also performed for SPC/E and TIP5P models of water. We prepared

model

〈τw〉 (fs)

στw (fs)

〈τj〉 (fs)

θj (deg)

∆F (kBT)

TIP3P SPC/E TIP5P

1304 1427 335

1767 1899 447

65 72 72

84 91 73

2.5 2.8 0.7

a 〈τw〉 is the mean waiting time, στw is the variance of the waiting-time distribution, 〈τj〉 is the jump time, θj is the magnitude of jump angle, and ∆F is the height of the free energy barrier.

similar systems of a (6,6) carbon nanotube solvated with these models of water. The confined water molecules were studied, and angular jumps were detected in all of them. The geometrical features (the changes in distances and angles) of the jump remained very similar for the various models of water studied. Figures 11 and 12 show the temporal evolution of the geometrical quantities, calculated across the jumps for TIP5P and SPC/E water under confinement. The results are summarized in Table 1. We observe that the three-point models TIP3P and SPC/E agree well in terms of the time scales associated with waiting and jumps, whereas TIP5P shows some difference. The mean waiting time for TIP5P is about one-fourth that of TIP3P. This is also reflected in the values of the free energy barrier, which are much larger for the TIP3P and SPC/E water models than for the TIP5P model. We attribute these differences to the fact that the charges for TIP5P are not

Reorientation of Water Molecules Confined in Nanotubes

Figure 13. Free energy barrier as a function of the relevant angular coordinate for the various models of water.

exactly located at the site of the oxygen atom, but rather are placed slightly behind it. This changes the potential energy landscape significantly, which leads to changes in the values of the time scales. The free energy landscape in the relevant angular coordinate is shown in Figure 13. 5. Conclusions Our MD simulations show that confined water molecules that are in a single-file arrangement inside narrow carbon nanotubes orientationally relax through angular jumps between two potential energy minima. This bistable behavior is very different from the mechanism of jumps in the bulk, where thermally activated jumps make the system explore many locally stable minima. The angular jump of a water molecule in the bulk involves the breaking of a hydrogen bond with one of its neighbors and the formation of a hydrogen bond with a different neighbor. In contrast, the angular jump of a confined water molecule corresponds to an interchange of the two hydrogen atoms that can form a hydrogen bond with the same neighbor. The free energy barrier between these two minima is a few kBT, with the electrostatic hydrogen-bonding energy accounting for most of this barrier. The waiting time between two jumps peaks at 60 fs, but it exhibits a long, slowly decaying tail. The value of the mean waiting time is ∼1 ps. The effect of jumps between these two minima shows up in the behavior of C1OH(t), which does not decay to 0 at long times. On the other hand, C2OH(t) does decay to 0 at long times. It should be mentioned that we considered here the orientational dynamics over time scales that are short compared to that of the reorientation of the dipole moment. Different behavior is expected for time scales over which the collective dipole moment flips. The analytic solution of a simplified model that qualitatively explains this “unusual” relaxation was also presented herein. Acknowledgment. We thank the Department of Science & Technology, India, for financial support. References and Notes (1) Ohmine, I.; Tanaka, H. Fluctuation, Relaxations, and Hydration in Liquid Water. Hydrogen-Bond Rearrangement Dynamics. Chem. ReV. 1993, 93, 2545–2566. (2) Marx, D.; Tuckerman, M. E.; Hutter, J.; Parrinello, M. The nature of hydrated excess proton in water. Nature 1999, 397, 601–604.

J. Phys. Chem. B, Vol. 113, No. 30, 2009 10329 (3) Bagchi, B. Water Dynamics in the Hydration Layer around Proteins and Micelles. Chem. ReV. 2005, 105, 3197–3219. (4) Nandi, N.; Bhattacharyya, K; Bagchi, B. Dielectric Relaxation and Solvation Dynamics of Water in Complex Chemical and Biological Systems. Chem. ReV. 2000, 100, 2013–2046. (5) Frank, H. S.; Wen, W. Y. Structure formation in liquid water. Discuss. Faraday Soc. 1957, 24, 133–140. (6) Stanley, H. E.; Teixeira, J. Interpretation of the unusual behaviour of H2O and D2O at low temperatures: Tests of a percolation model. J. Chem. Phys. 1980, 73, 3404–3422. (7) Belch, A.; Rice, S. A. The distribution of rings of hydrogen-bonded molecules in a model of liquid water. J. Chem. Phys. 1987, 86, 5676– 5682. (8) Rahman, A.; Stillinger, F. H. Molecular Dynamics Study of Liquid Water. J. Chem. Phys. 1971, 55, 3336–3359. (9) Luzar, A.; Chandler, D. Hydrogen bond kinetics in liquid water. Nature 1996, 379, 55–57. (10) Luzar, A. Resolving the hydrogen bond dynamics conundrum. J. Chem. Phys. 2000, 113, 10663–10675. (11) Laage, D.; Hynes, J. T. A Molecular Jump Mechanism of Water Reorientation. Science 2006, 311, 832–835. (12) Laage, D.; Hynes, J. T. Reorientational Dynamics of water molecules in anionic hydration shells. Proc. Natl. Acad. Sci. 2007, 104, 11167–11172. (13) Laage, D.; Hynes, J. T. Do more strongly hydrogen-bonded water molecules reorient more slowly. Chem. Phys. Lett. 2006, 433, 80–85. (14) Ivanov, E. N. Theory of Rotational Brownian Motion. SoV. Phys. JETP 1964, 18, 1041. (15) Stenger, J.; Madsen, D.; Hamm, P.; Nibbering, E. T. J.; Elsaesser, T. Ultrafast Vibrational Dephasing of Liquid Water. Phys. ReV. Lett. 2001, 87, 027401-1–027401-4. (16) Woutersen, S.; Emmerichs, U.; Bakker, H. J. Femtosecond MidIR Pump-Probe Spectroscopy of Liquid Water: Evidence for a TwoComponent Structure. Science 1997, 278, 658–660. (17) Lawrence, C. P.; Skinner, J. L. Vibrational spectroscopy of HOD in liquid D2O. III. Spectral diffusion, and hydrogen-bonding and rotational dynamics. J. Chem. Phys. 2003, 118, 264–272. (18) Woutersen, S.; Bakker, H. J. Ultrafast Vibrational and Structural Dynamics of the Proton in Liquid Water. Phys. ReV. Lett. 2006, 96, 138305(1)138305(4) . (19) Fecko, C. J.; Eaves, J. D.; Loparo, J. J.; Tokmakoff, A. A.; Geissler, P. L. Ultrafast Hydrogen-Bond Dynamics in the Infrared Spectroscopy of Water. Science 2003, 301, 1698–1702. (20) Zangi, R. Water confined to a slab geometry: A review of recent computer simulation studies. J. Phys.: Condens. Matter 2004, 16, S5371– S5388. (21) Tan, H. S.; Piletic, I. R.; Fayer, M. D. Orientational dynamics of water confined on a nanometer length scale reverse micelles. J. Chem. Phys. 2005, 122, 174501 (1-9). (22) Zhu, F.; Schulten, K. Water and Proton Conduction through Carbon Nanotubes as Models for Biological Channels. Biophys. J. 2003, 85, 236– 244. (23) Hanasaki, I.; Nakamura, A.; Yonebayashi, T.; Kawano, S. Structure and stability of water chain in a carbon nanotube. J. Phys.: Condens. Matter 2006, 20, 015213-1–015213-7. (24) Best, R. B.; Hummer, G. Reaction Coordinates and Rates from Transition Paths. Proc. Natl. Acad. Sci. U. S. A. 2005, 102, 6732–6737. (25) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. Comparison of simple potential functions for simulating liquid water. J. Chem. Phys. 1983, 79, 926–935. (26) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. The missing term in effective pair potentials. J. Phys. Chem. 1987, 91, 6269–6271. (27) Mahoney, M. W.; Jorgensen, W. L. A five-site model for liquid water and the reproduction of the density anomaly by rigid, nonpolarizable potential functions. J. Chem. Phys. 2000, 112, 8910–8922. (28) Case, D. A.; Pearlman, D. A.; Caldwell, J. W.; Cheatham T. E., III; Wang, J.; Ross, W. S.; Simmerling, C. L.; Darden, T. A.; Merz, K. M.; Stanton, W. A.; Cheng, A.; Vincent, J. J.; Crowley, M.; Tsui, V.; Gohlke, H.; Radmer, R.; Duan, Y.; Pitera, J.; Massova, I.; Seibel, G. L.; Singh, U. C.; Weiner, P.; Kollman, P. A. AMBER 7; University of California: San Francisco, CA, 1999. (29) Mukherjee, B.; Maiti, P. K.; Dasgupta, C.; Sood, A. K. Strong correlations and Fickian water diffusion in narrow carbon nanotubes. J. Chem. Phys. 2007, 126, 124704-1–124704-8. (30) Hummer, G.; Rasaiah, J. C.; Noworyta, J. P. Water conduction through the hydrophobic channel of a carbon nanotube. Nature 2001, 414, 188–190. (31) Laage, D.; Hynes, J. T. On the Molecular Mechanism of Water Reorientation. J. Phys. Chem. B 2008, 112, 14230–14242. (32) James, T.; Wales, D. J. Energy Landscapes for Water Clusters in Uniform Electric Field. J. Chem. Phys. 2007, 126, 054506-1–054506-13.

10330

J. Phys. Chem. B, Vol. 113, No. 30, 2009

(33) Berne, B. J.; Pecora, R. Dynamic Light Scattering with Applications to Chemistry, Biology and Physics; Wiley: New York, 1976. (34) Winkler, K.; Lindner, J.; Bursing, H.; Vohringer, J. Ultrafast Ramaninduced Kerr effect of water: Single molecule versus collective motions. J. Chem. Phys. 2000, 113, 4674–4682. (35) Lombardo, T. G.; Debenedetti, P. G.; Stillinger, F. H. Computational probes of molecular motion in the Lewis and Wahnstrom model for orthoterphenyl. J. Chem. Phys. 2006, 125, 174507. (36) Seki, K.; Bagchi, B.; Tachiya, M. Orientational relaxation in a dispersive dynamic medium: Generalization of the Kubo-Ivanov-Anderson jump-diffusion model to include fractional environmental dynamics. Phys. ReV. E 2008, 77, 031505-1–031505-10.

Mukherjee et al. (37) Mukherjee, B.; Maiti, P. K.; Dasgupta, C.; Sood, A. K. Strongly anisotropic orientational relaxation of water molecules in narrow carbon nanotubes and nanorings. ACS Nano 2008, 2, 1189–1196. (38) Loughnane, B. J.; Scodinu, A.; Fourkas, J. T. Extremely Slow Dynamics of a Weakly Wetting Liquid at a Solid/Liquid Interface: CS2 Confined in Nanoporous Glasses. J. Phys. Chem. B 1999, 103, 6061– 6068. (39) Loughnane, B. J.; Fourkas, J. T. Geometric Effects in the Dynamics of a Nonwetting Liquid in Microconfinement: An Optical Kerr Effect Study of Methyl Iodide in Nanoporous Glasses. J. Phys. Chem. B 1998, 102, 10288–10294.

JP904099F