Open Access is an initiative that aims to make scientific research freely available to all. To date our community has made over 100 million downloads. It’s based on principles of collaboration, unobstructed discovery, and, most importantly, scientific progression. As PhD students, we found it difficult to access the research we needed, so we decided to create a new Open Access publisher that levels the playing field for scientists across the world. How? By making research easy to access, and puts the academic needs of the researchers before the business interests of publishers.

We are a community of more than 103,000 authors and editors from 3,291 institutions spanning 160 countries, including Nobel Prize winners and some of the world’s most-cited researchers. Publishing on IntechOpen allows authors to earn citations and find new collaborators, meaning more people see your work not only from your own field of study, but from other related fields too.

The SEC Strategic Research Cluster and the Centre for Synthesis and Chemical Biology, School of Chemical and Bioprocess Engineering, University College Dublin, Belfield, Dublin 4, Ireland

José-Antonio Garate

The SEC Strategic Research Cluster and the Centre for Synthesis and Chemical Biology, School of Chemical and Bioprocess Engineering, University College Dublin, Belfield, Dublin 4, Ireland

J. M. Don MacElroy

The SEC Strategic Research Cluster and the Centre for Synthesis and Chemical Biology, School of Chemical and Bioprocess Engineering, University College Dublin, Belfield, Dublin 4, Ireland

*Address all correspondence to:

1. Introduction

The accelerated appearance of applications and theoretical studies of carbon nanotubes (CNTs)^{1} in the literature over the last two decades reveals their broad attraction in fields as diverse as nanoelectronics, bionanotechnology and solar energy conversion nanodevices^{2-5}. In particular, they provide the basis for simplified models of nanofluidic machines and membrane channels such as aquaporins^{6-9} and have been employed in a range of simulation studies at the nanoscale to investigate the transport of water undergoing single-file or highly confined diffusion and permeation in a single CNT or in arrays of CNTs^{10-22}.

The simulation method of choice in the study of water transport through nanostructures is molecular dynamics (MD) (see, for example, [23-25]) which provides real-time trajectories of all of the atoms within the system under investigation and, through the statistical framework of linear response theory, permits a direct link with macroscopic observables. One of the first studies of water confined within CNTs employing MD, supplemented with an analysis using a continuous-time random walk model (CTRWM), was reported by Berezhkovskii and Hummer^{12,13} and demonstrated that water fills a (6,6) CNT spontaneously in a single-file fashion. In later work, MD simulations focused on pressure-driven osmotic water flow in single CNTs and arrays of CNTs^{14-16} and in model pores containing charge distributions which are known to be ubiquitous for biological channels^{17,18}. Water channels made of double-walled carbon nanotubes embedded in a lipid membrane, where the inner tube is modeled as a (6,6) CNT with a charge distribution emulating the NPA motif of aquaporins^{19}, have also been investigated, and anomalous observations relating to the dynamics, configurational states and apparent pore-length dependence of the properties of water confined within CNTs have received attention^{20-22}.

In [26,27], molecular simulations were reported for unidirectional single-file diffusion of water within (5,5) and (6,6) single-walled carbon nanotubes (SWCNT) and normal diffusion in (8,8) and (11,11) SWCNT’s across a lipid membrane. These studies employed both equilibrium MD simulation techniques and an adaptation of MD to investigate conditions involving the influence of both static and alternating electric fields. While the fields employed in this work were relatively weak vis-à-vis fields generated by charge distributions, it was observed that the permeation of water for the (5,5) SWCNT embedded in a lipid membrane is enhanced by static electric fields due to a decrease in the fluctuations of the number of water molecules inside the SWCNT. In the larger single-file (6,6) CNT pore, dipole alignment with the imposed external field is more easily achieved due to a greater degree of rotational freedom of the water molecules and this was observed to result in a significant reduction in the water self-diffusion flux^{26}. In this work it was also observed that the non-field (6,6) flux results were much higher than those reported by Zhu and Schulten^{14}, even though the same potential/force field was used in both studies. We suggested that this difference may be due to the array of parallel nanotubes employed by Zhu and Schulten which could cause friction (via long-range electrostatics) between the single files of water. This was confirmed in [27] where it was demonstrated that in these CNT-arrays there is a significant friction effect that retards movement of the single files of water.

In this chapter, we report further studies of non-equilibrium MD simulations of water permeation within (6,6) SWCNT’s across a lipid bilayer under the influence of both static and alternating electric fields, taking advantage of GPU-accelerated MD^{28} to extend the simulation timescales to longer than was previously possible^{26,27}, while also investigating single-molecule translational and rotational frequency modes, collective properties and friction in tube-arrays in greater detail, allied with straightforward, mechanistic models to interpret our results.

All of the GPU-accelerated MD simulations^{23-25} were carried out using the NAMDv2.7^{28,29} program suite, in conjunction with the CHARMM27 parameter set^{30} and the TIP3P model for water^{31}. Initially, a periodically replicated single-nanotube system was studied, comprising a (6, 6) armchair SWCNT, with a diameter of 8.20 Å and approximate length of 36.9 Å, embedded in a 1-palmytoil-2-oleoyl-sn-glycero-3-phosphatidylcholine (POPC) membrane patch. Layers of water of around 10Å were added in the +z and -z-direction to solvate the system. Each atom of the SWCNT was modelled as a neutral sp^{2} aromatic carbon with uncapped ends (in a similar way as shown in Fig. 1, although this is for arrays of four CNTs). For all of these systems, the plane of the POPC membrane was placed in the x-y plane, and the axis normal to the membrane was along the z-direction, with the CNT oriented along this z-axis (perpendicular to the plane of the membrane – cf. Fig. 1). The final size of the (three-dimensional-periodic) simulation box was 34.84 x 34 x 66.29 Å^{3}, with total atom numbers of 8,415 (1337 water molecules and 30 phospholipids).

To model the arrays of CNTs, three-dimensional-periodic cells consisting of four (6, 6) SWCNT arrays were constructed with separations of 15, 20 and 25 Å between the geometric central axes of each (see Fig. 1), and inserted along the z-axis of a POPC membrane, again solvated with layers of water molecules of around 10 Å in thickness in the +z and -z-directions. In each of these cases, the membrane was also placed in the x-y plane, with the principal axis of the SWCNT oriented perpendicular to the membrane. The final sizes for the simulation boxes in the array-systems were 45.75 x 60.32 x 66.68 Å^{3}, 49.36 x 49.21 x 71.27 Å^{3} and 63.10 x 59.99 x 66.27 Å^{3} for respective separations of 15, 20 and 25 Å, with total atom numbers of 18,979 (3045 water molecules and 62 phospholipids), 17,901 (3043 water molecules and 54 phospholipids) and 26,309 (4327 water molecules and 88 phospholipids), respectively. All molecular systems were constructed using the VMDv1.86^{32} program.

The Particle Mesh Ewald^{33} method was utilised to treat long-range electrostatics to within a relative tolerance of 1x10^{-6} with a cut-off distance of 12 Å applied to real-space Ewald interactions; it was found that a real-space screening coefficient of 0.258 Å^{-1} and (three-dimensional Fourier transform) grid spacings in the 0.9 to 1 Å range in each Cartesian direction led to this desired accuracy level, along with close-to-optimal execution speeds. A cut-off distance of 12 Å was also applied to van der Waals interactions, using a smooth switching function between 10 and 12 Å. Using the r-RESPA method^{34}, multiple time steps were applied to MD, with 1 fs for bonded interactions, 1 fs for short-range non-bonded interactions (including real-space Ewald interactions) and 2 fs for the full (real- and reciprocal-space) electrostatics evaluation. All production runs were carried out using Langevin dynamics^{35} coupled to an NVT reservoir set at 300 K with a damping coefficient of 1 ps^{-1}, representing rather mild coupling which affects diffusion little. The SHAKE^{36} algorithm was used to constrain bond lengths to all hydrogen atoms.

Prior to production, an NPT relaxation run was carried out for 1 ns; a Nosé-Hoover approach was applied^{37} using anisotropic cell variation, with Langevin dynamics for piston fluctuation control^{38}; here, the set points were 1 atm and 300 K. During this NPT relaxation, the nanotube atoms were fixed, while the lipid bilayer tails and water molecules were allowed to move freely. Following this, the CNT and lipid-head atoms were released slowly, keeping a harmonic constraint thereon with a strength of 1 kcal/mol Å². This strength coefficient was reduced progressively to 1.5, 1.0 and 0.5 kcal/mol Å² every 100 ps for 0.4 ns. After this initial relaxation, during which the pressure and volume had stabilised, MD simulations of each of the (6, 6) CNT-systems (single-tube and array) were carried out for 50 ns in the NVT ensemble (in the absence and presence of external fields – vide infra). Three replicates were performed for each system in each case (zero- and in-field).

Use of GPU-accelerated MD implemented in NAMD v2.7^{28,29} on an NVIDIA double-precision C2050 platform led to a substantial acceleration (approximately 2.5-fold) vis-à-vis the equivalent number of CPU-cores on quad-core Intel Xeon Harpertown-Stoakley nodes connected via sub-microsecond-latency Infiniband. The acceleration was primarily, though not only, in the evaluation of intermediate-range non-bonded van der Waals and electrostatics interactions, where the computational overhead was reduced approximately six-fold relative to that required with the equivalent number of CPU-cores in a CPU-only approach (on the abovementioned Intel Xeon chips).

Static electric fields (EFs) were applied in the non-equilibrium MD simulations with an intensity of 0.0065 V/Å in the +z-direction (i.e., along the CNT-axis). The EF exerts a force over partial atomic charges ia:

fia= qiaEE1

Non-equilibrium MD simulations were also carried out in the presence of alternating electric fields^{39} along the same z-axis (i.e., along the axes of the nanotubes). In this case, the RMS intensity was set to 0.0065 V/Å, the same as that of the peak frequency for the static EF, but the field does change direction every period between the +z and –z axis, ipso facto^{40,41}. Three different frequencies were used: 2.45, 50 and 100 GHz. Here, the force exerted over the partial atomic charges is given by:

fia=qiaE0cos(ωt)E2

No constraints were applied to the membrane during MD in either the absence or presence of static and alternating electric fields, while a light 0.5 kcal/mol Å^{2} harmonic constraint was applied to the SWCNT atoms.

3. Modelling diffusion and permeation of water in single-file channels

The transport of water through single file channels embedded in lipid bilayers may be analysed from two perspectives (a) self-diffusion of water (e.g. labeled isotopically for direct physical measurement or, for simulation purposes, simply assigned a ‘colour’) and (b) water permeation under the influence of a pressure gradient (osmotic or otherwise). The two parameters which characterise these transport processes are the diffusive and osmotic permeabilities, p_{
d
} and p_{
f
}, respectively (see [10] for a review of details). These are deﬁned through the flux expressions:

jd=pdΔctrE3

and

jf=1kBTpfΔPE4

where Δctr is the tracer (or colour) concentration difference across the membrane, P is the trans-membrane liquid pressure difference (due, under normal circumstances, to an osmotic bulk solute concentration difference), T is absolute temperature and k_{
B
} is Boltzmann’s constant.

In [10] and references cited therein (see also [12-15]) the diffusive permeability has been investigated using a continuous time random walk model (CTRWM) and, in the absence of a chemical potential difference Δμ or an external force, the osmotic permeability was computed by equilibrium MD using a collective diffusion model. In [26, 27], a more general model was developed for the diffusive permeability (which is not restricted in application to single-file pores) and this more general approach is described below. In addition, the collective transport permeability is analysed in this chapter in greater detail using the generalized linear transport theory for slip flow in pores developed in earlier work^{42,43}.

3.1. The diffusive permeability

To compute the diffusive permeability, p_{
d
}, of water in a strongly inhomogeneous system constituting two bulk liquid reservoirs separated by a lipid membrane with embedded carbon nanotubes, the simplest and statistically most accurate approach is to employ a colour diffusion technique. In this approach the flux of ‘black’ water molecules diffusing from left to right in Fig. 2 is monitored during the course of a long molecular dynamics simulation (a similar analysis applies to ‘white’ water molecules diffusing in the opposite direction). If the cross-sectional areas in the bulk and in the pore are A and A_{
p
}, respectively, and L is the length of the mixing zone between the black reservoir on the left of Fig. 2 and the pore mouth, then at steady state the net flow (particles per unit time) between the black reservoir and the outer perimeter of the CNT is given by

jd=DAL(cB(1)−cB(2))E5

where c_{
B
}(1) is the concentration of black water molecules in the black reservoir and c_{
B
}(2) is the concentration of black water molecules just outside the CNT pore mouth.

At the pore mouth a simple reversible rate process takes place for the molecules to either enter the pore or return into the reservoir on the left and the net rate is given by

r2−2p=kfcB(2)−krcBp(2)E6

where the forward and reverse rate constants k_{
f
} and k_{
r
} are in units of sec^{-1} and c_{
Bp
}(2) is the concentration of black water molecules just inside the CNT pore adjacent to the pore mouth. The net rate in Eqn. (6) is a volumetric measure (these rate constants will be determined to a great extent by the hydrogen bond structure and dynamics). Multiplying by a volume V_{i} (an interfacial volume corresponding to A_{
p
}δ with δ = length scale of the order of the molecular dimensions representing the scale of the concentration inhomogenities) then for steady state conditions we have

where K_{
p
} is defined as k_{
f
}/k_{
r
} and p_{1} is the partial permeability for diffusive transport from the black reservoir to just inside the CNT:

1p1=LDA+1kfViE9

The steady-state flow along the pore may also be defined using an intra-pore diffusion coefficient, D_{
p
} (which is considered to exist in the same sense as proposed in [10, 12-14] but which has a deeper meaning than indicated by the simple CTRWM)

jd=DpApLp(cBp(2)−cBp(3))E10

where L_{
p
} is the length of the pore and c_{
Bp
}(3) is the concentration of black water molecules just inside the CNT pore adjacent to the pore mouth near the white particle reservoir.

Also, for the right hand side of the system it follows from Eqn. 7 that

jd=p1Kp(cBp(3)−KpcBp(4))E11

where c_{
Bp
}(4) is the concentration of black water molecules in the white particle reservoir.

Now combining Eqns. 8, 9 and 10 eliminating the intermediate concentration terms gives

jd=pd(cB(1)−cB(4))E12

where

1pd=2LDA+2kfVi+LpDpKpApE13

Eqn. 8 simply states that the permeability is governed by a sum of mass transfer resistance terms due to bulk diffusion, incorporation into the pore mouth and intra-pore diffusion. Normally, one sets c_{
B
}(4) = 0 (in the simulations the ‘black’ water molecules entering the ‘white’ reservoir are simply switched in colour to ‘white’) and c_{
B
}(1) = molecular/molar density of water in the bulk i.e. c_{
H2O
} = 33.5 nm^{-3} or 55.6 mol.litre^{-1}. Therefore

jd=pdcH2OE14

The overall permeability can be simply determined from the MD simulations and by selecting different positions for the ‘black’ and ‘white’ reservoirs, one can eliminate the bulk diffusion D/L terms.

In the development of the CTRW model reported by Berezhkovskii and Hummer^{12}, there is the implicit assumption that L is sufficiently small to minimise bulk-phase diffusional resistance while retaining a sub-nanometre range close to the pore mouth (δ ~ a) to ensure that entry/exit kinetics are fully taken into consideration. This condition with the additional assumption that the single file is considered to be densely packed leads to

Lp=Na;δ=aE15

where a is a length scale characteristic of the size of the water molecules. Furthermore, the diffusion coefficient within the pore is given by the CTRW model as

Dp=a22τ=a22krE16

In our model, we employ the physically reasonable assumption that the intra-pore hopping rate k is equal to k_{
r
}, i.e., the reverse rate constant in Eqn. 5. Within the pore, the forward and reverse hopping rate constants are equal. At the pore mouth, however, we distinguish between the forward (pore entry) and the reverse (pore exit) rate constants. For the CNT’s under consideration here, one may anticipate k_{
r
}> k_{
f
} or K_{
p
}< 1 and if we assume the molecular volumevw=aAp, then we find

Permeation rate=kfaAp21N+1cH2O=vwkoN+1cH2OE17

which is the result provided by the CTRW model in [12]. This analysis demonstrates the level of the assumptions involved in the CTRWM approach. As noted above, in a more general setting, the in-pore forward hopping frequency is not the same as at the pore mouth. Furthermore, the more general model of Eqn. 11 applies to wide (i.e. multi-pass) pore conditions of arbitrary cross-sectional area and it is also simple to extend Eqn. 11 to include rate processes at other junctions (e.g. local charged sites) within the pore and to relate the rate constants to the energy barriers associated with the particle movement along the channel.

3.2. The osmotic permeability, collective properties and friction in pore-arrays

In general, the flux of water molecules across a membrane subject to a chemical potential driving force, Δμ, which in turn originates from an osmotic pressure gradient, may be expressed as

jf=cH2OkBTpfΔμE18

which is equivalent to Eqn. (4). Here, cH2Ois the molecular density of water while p_{
f
} is osmotic permeability as defined, inter alia, in [10, 14-16]. Using results from linear response theory^{42,43} the permeability of membranes investigated in this work may be expressed as

wherecH2O,M is the water concentration within the membrane (which is related to cH2O via the partition coefficient^{44}K between the bulk liquid phase on both sides of the membrane and the membrane pores), N_{
T
} is the total number of water molecules within the membrane at any given time and DM(z) is the axial slip-diffusion (or viscous slip) coefficient of the water molecules within the CNT-pores of the membrane. The function in the integral of Eqn. 12 is the autocorrelation function (ACF) of the axial, collective centre-of-mass velocity of all of the water molecules present in the pores of the membrane, with the collective (number averaged) velocity given by

uz(t)=1NT∑i=1NTvi,z(t)E20

where vi,z(t) is the velocity of water molecule i along the z-axis of the CNT-pores. In Eqn. 12, z(t) is the corresponding collective, axial centre-of-mass position of the water molecules in the pore at any point in time:

z(t)=1NT∑i=1NTzi(t)E21

For a membrane containing an array of N_{
p
} pores, each of which is accessible to water to the same extent, as studied in this work, we note that Eqn. 12 may be expressed as:

where N(j) is the number of water molecules in pore j. The correlation function in Eqn. 16 may be considered a total one (TCF), in that it applies to all water molecules in the membrane (i.e. present in all of the CNT-pores). Considering Eqn. 16, it is clear that this may be written in terms of the ACFs of the collective velocities in individual pores as well as cross-correlation functions (CCFs) of collective velocities between particular pores. Eqn. 16 demonstrates the possibility of significant cross-interactions between the individual pores in the array, unless they are well separated. Should this latter condition apply, then Eqn. 16 simplifies to:

where the subscript ∞ denotes the limit of infinite separation of all pores, and pf* is the permeability for a single, isolated pore.

In order to explore further the underlying nature of the collective pore-velocity correlations of Eqn. 16 which govern the distance-dependent friction observed in our previous work^{27} it is appropriate to consider a simplified model of friction. For a given collective velocity within the pores, uz, the friction may be approximated by a term due to isolated pores (the ‘∞’ case), as well as an additional, distance-dependent term (a finite separation ‘f’ case), which will be zero in the case of the limit of infinite separation:

duzdt=−(1/τ∞+1/τf)uzE24

the solution of which is

uz(t)=uz(0)exp(−(1/τ∞+1/τf)t)E25

Here, the isolated-case 1/τ∞ friction term may be considered as a relaxation time which might possibly be expected in the decay of the ACF of the isolated-case collective velocity, while that of 1/τf would be additional friction that might be manifest in the decay of the TCF of Eqn. 16, vis-à-vis that of the isolated case. At its simplest, one might model the distance-dependent friction term, 1/τf, to be inversely proportional to separation d, i.e.1/τf~ 1/d. Exploring this collective-velocity correlation further for the CNT-array:

4.1. Static and alternating electric field effects and water self-diffusion

The permeation rates for water molecules through the SWCNT (in molecules per ns) embedded within the lipid layer at steady-state, in the absence and presence of axially-applied static and alternating electric fields, are provided in Fig. 3. Although the rate in the static field is marginally higher (notwithstanding overlapping error bars vis-à-vis the zero-field case), the alternating field results show no discernible trend when the electric field frequencies are incremented. It is important to note that each point represents the average fluxes of three independent simulations. Statistically, with the possible exception of the static-field case, we observe that there is no correlation between the field frequencies and the computed water permeation rates. This is in essential agreement with our previous study^{27} although longer MD timescales from GPU-acceleration has led to reduced error bars in this work. In our previous work^{26,27} monotonic decreases of water permeation rates with field frequencies were found for axial field application to (5,5) CNT-systems towards the zero-field rate at higher frequency, but not to (6,6) CNT-cases (or larger CNTs). For (5,5) CNT-systems, it was found that axial static electric fields enhanced water self-diffusion across them due to a decrease in the fluctuations of the number of water molecules within the pores, arising from the water molecules’ dipole orientation with the field^{26}. Although it is difficult to reach a conclusion for the (more statistically accurate vis-à-vis [27]) slightly larger static-field rate in the case of this work for (6,6) CNTs, due to overlapping error bars with respect to the zero-field result, the single-file nature of water transport in (6,6) CNTs, similar to the (5,5) case, suggests that a more minor decrease in the fluctuations of the number of water molecules inside is at play in the static field. This was observed in this study with results slightly lower than previous (6,6) data^{26,27} although error bars for the size of the fluctuations were again overlapping vis-à-vis the zero-field case, so no firm statistical conclusion can be drawn in that regard. It would be necessary to carry out a greater number of longer runs to investigate this further, somewhat prohibitive even with ample GPU-accelerated computing resources.

The water dipole orientations within the pore (see Fig. 4), illustrate more clearly the influence of the external field on the probability distributions of the cosine of the angle of the dipole with the positive z-axis, α. In agreement with our earlier results, the statistically more accurate results reported here illustrate that at high frequencies the distributions become more symmetric with peaks at -1 and 1, in very good agreement to the zero-field case. In the static field case the distribution is shifted almost entirely to the left side of the diagram.

From consideration of the fluxes and dipolar orientations, it is clearly evident that the results of the low-frequency electric field simulations are more similar to the static-field case, whereas high-frequency fields resemble zero field conditions. These results support the finding that the effects of low-frequency electric fields should be closer to those of static electric fields (if any), because the slower time variation of the field (i.e., longer period) relative to permeation and loading timescales allows the field to influence the dynamics of these events. On the other hand, high-frequency electric fields do not readily permit enhancements in permeation or alterations in dipole alignment vis-à-vis the zero-field case, due to their faster oscillation and shorter period relative to permeation timescales in (6,6) CNTs.

Considering single-molecule translational and librational modes in the isolated CNT system, the power spectra of the oxygen- and hydrogen-atom VACFs are shown in Figs. 5 and 6, respectively, under zero-field conditions, along with corresponding results for bulk water. The spectra are decomposed in terms of x-y and z-directions, and a translational ‘bumping’ mode vis-à-vis the CNT walls is clearly visible in the x-y plane at around 35 cm^{-1} (circa 0.8-0.9 ps) in the oxygen-atom spectra, along with corresponding librational modes of 35 and 160 cm^{-1} (around 0.18 and 0.8-0.9 ps). The axial-spectra, although somewhat more bulk-like, also exhibit important differences vis-à-vis the bulk state: notably a greater translational ‘rattling’ at 260 cm^{-1} (around 0.12 ps), and a marked reduction in higher-frequency librational motion, in accord with retarded rotational motion (evident in the preferred dipole distributions of Fig. 4).

4.2. Water permeation through nanotube-arrays

Steady-state water diffusion rates, j_{
d
}, and hopping rates, k, were computed for the (6,6) SWCNT array simulations and are reported in Table I. The intra-pore water self-diffusion coefficients, D_{
p,} provided in this Table were computed independently as described in [26] using the Einstein relation for self-diffusion (a special case of Eqn 12)

Dp=limt→∞12d(Δz)2dtE29

The results summarized here are essentially in agreement with previously reported values^{25,26} albeit with reduced error bars due to the longer timescales rendered possible from GPU-acceleration. There is an incremental increase in j_{
d
} when the separation between the pores is augmented, which is a direct consequence of the incremental increase in D_{
p
} with the separation of the nanotubes for separations greater than 15 Å. The lower than anticipated diffusive flux for the 15 Å close packed structure when comparing the magnitude of j_{
d
} relative to D_{
p
} arises in part due to the distance-dependent friction effect^{27} that retards water movement within the pores as demonstrated by the osmotic permeability results (the right-hand column of Table I) but also reflects a greater degree of interference between the pores at the pore mouths, namely a lowering of the forward (pore entry) rate constant k_{
f
} in the second term in Eqn. 11.

Estimates of the hopping rate constant, k, suggested by the CTRW model were also computed to compare our results with those reported in refs. [12-14], since, as previously shown^{26,27} our results for the zero-field (6,6) CNT were in quite good agreement with those reported by Berezhkovskii and Hummer^{12,13} but differed with the ones reported by Zhu and Schulten^{14}, even though both our results and the ones previously mentioned are based on the TIP3P/CHARMM potential. The hopping rate result, k = 57.4 ns^{-1}, reported here for simulations involving the 15 Å separation is more than a factor of 2 higher than the corresponding value reported by Zhu and Shulten (k = 26.9 ns^{-1}). As discussed in [27] only four tubes were employed in our simulations (embedded in a non-diffusing phospholipid membrane) while a close packed periodically imaged hexagonal array of twelve tubes was employed in the work conducted by Zhu and Schulten.^{14} The much lower k observed in ref. 14 is explained by an increase in the friction effect generated by the twelve single files of water interacting via long range electrostatics, rather than the four nanotubes employed in our work.

Separation (Å)

j_{d}(ns^{-1})

p_{d}(nm^{3}.ns^{-1})

D_{p }(nm^{2} ns^{-1})

^{a}k (ns^{-1})

^{b}p_{f} (nm^{3}.ns^{-1})

^{c
p
f
F
M
}(nm^{3}.ns^{-1})

15

1.69 (0.18)

0.050

1.94 (0.07)

57.40

0.896

1.74

20

2.12 (0.28)

0.063

2.00 (0.06)

59.17

1.156

1.85

25

2.32 (0.13)

0.069

2.14 (0.11)

63.31

1.268

1.93

^{
∞
}

2.47 (0.29)

0.074

2.38 (0.12)

70.41

1.368

2.18

Table 1.

The hopping rates defined by the CTRW model, k (≡ k_{
r
}) = 2D_{
p
}/a^{
2
}, where a = 0.26 nm is the projected length for a given water molecule along the pore axis (for more details see [12]).^{b} The osmotic permeability predicted by integration of the viscous slip total velocity-CF (TCF) (Eqn. 16).^{c} The osmotic permeability predicted by the simplified friction model (Eqn. 19).Steady-state tracer permeation rates, tracer intrapore diffusion coefficients, hopping rates and osmotic permeabilities for the CNT arrays.

To demonstrate the existence of significant frictional effects between the aligned water files in [27] we estimated an effective viscosity along the axis of the water confined inside the CNTs. In summary, this was done by use of the appropriate Einstein relation^{23} applied to tracking chains of six water molecules initially nearest the centre of each of the four nanotubes. The four hexamers of water so defined were treated as a fluid of four ‘particles’ and their motion was monitored until one or more of the water molecules in any given hexamer left its respective nanotube. The hexamers were then redefined and sampling continued until the next disruption in one of the hexamer chains. Averaging over each such occurrence in all four tubes in the zero-field simulations then enabled the determination of the effective viscosity of the files using:

η=V2kBTlimt→∞ddt(Pzz(t)−Pzz(0))2 where Pzz=1V∑i=1NTrizpizE30

The long-time slopes of the mean square displacement of the dynamical variable Pzz as functions of time (see Fig. 7) for each nanotube separation provided the effective viscosities 7.5 ± 0.45 (15 Å), 5.5 ± 0.4 (20 Å) and 4.4 ± 0.32 (25 Å) × 10^{-3} Nsm^{-2}. These are agreement with our previous study^{27} although the standard deviations are reduced due to longer simulations. The significantly larger effective viscosity the closer the tube-tube separation confirms the strong interaction induced by the surroundings on the water molecules confined within each individual CNT pore. The dependence of η on the separation distance between the nanotubes is more clearly demonstrated in Fig. 8 where the fit to a simple 1/d dependence is shown to be excellent.

To investigate friction further, collective, axial VACFs of all of the water molecules in each individual CNT along the z-axis are shown in Fig. 9a, including the isolated tube case. A log-plot of this in Fig. 9b shows an essentially negative-exponential decay, along the general trend of an Aexp(−t/τ)cosωt motif. Fitting to the successive peaks in the log-plot yields an isolated tube frictional relaxation time of 0.16 ps. From the other slopes, the distance-dependent frictional contribution may be determined for each of the four-CNT arrays (inferred from the fitted peak lines), and these are detailed in Table II. The 1/τ_{f} distance-dependent friction parameter declines with increasing separation, and a plot of this versus the inverse separation as shown in Fig. 8, coupled with a regression line, shows an approximate 1/d relationship. Direct evaluation of the z-axis collective-velocity cross correlation function (VCCF) for each tube with the others in the three CNT-array cases confirms that the magnitude of the VCCF declines with increasing separation, and this is shown in Fig. 10.

The osmotic permeabilities were evaluated in two independent ways, and the results are provided in the last two columns of Table II. Integration of the Green-Kubo-type representation of the viscous slip total velocity correlation function (TCF) (Eqn. 16) (see [42,43] and also [23, 45,46] for additional details on Green-Kubo analysis) leads to convergence within 10-15 ps and the results for the osmotic permeability are provided in the second-last column of Table II. These results are the exact values for the osmotic permeabilities. The larger values obtained from the approximate friction analysis (Eqn. 19) are to be expected given that the friction model does not take into account negative backscattering parts of the velocity correlation function clearly observed in Fig. 9a. Nonetheless, the general trend with nanotube separation is clearly evident. Also note that the ratio of the osmotic and diffusion permeabilities p_{
f
}/p_{
d
} in all cases (with p_{
f
} given by the exact results in the second last column of Table I) is approximately 18 in magnitude. This corresponds reasonably with the estimate of 15 obtained from an approximation suggested in reference [10].

Separation (Å)

τ_{eff} (ps)

τ_{f} (ps)^{*}

15

0.128

0.64

20

0.134

0.825

25

0.142

1.26

^{
∞
}

0.16

Infinite

Table 2.

Inferred from 1/τ_{eff} = 1/τ_{eff,}∞ + 1/τ_{f} Table II. Fitted values for effective frictional relaxation times, τ_{
eff
}, and values of the resultant separation-friction relaxation time, τ_{
f
}. See plot of 1/τ_{
f
} versus 1/d (Fig. 8).

Single-file diffusion results through (6,6) SWCNTs, embedded in a lipid bilayer, under the influence of static and time-varying electric fields, have been reported. No substantial influence of alternating electric fields on water flux through SWCNTs of larger diameters was observed, although dipole orientations were affected markedly, with higher-frequency distributions resembling those of the zero-field case. In the case of static fields, however, the water flux appears somewhat higher than that of the zero-field case, owing possibly to a similar reduction in water-number fluctuations therein, as found in previous work for (5,5) CNTs^{26,27}; an overlapping of error bars vis-à-vis the zero-field case, however, renders any firm statistical conclusion elusive in this respect. The (RMS) field intensities of 0.0065 V/Å are weak compared to previous molecular simulation studies employing static or time-varying electric fields^{39} and also weak relative to fields generated by charge distributions^{14}, and are substantially more realistic with respect to experiment.

When arrays of CNTs are employed in the simulations the interpretation on the transport properties of the fluid within the individual nanotubes must take cognisance of inter-tube interactions. The (6,6) CNT-array simulations reported in this work have demonstrated the presence of a friction effect due to long-range electrostatic interactions between the moving water single-files within the pores, with better statistics than was possible in previous work.^{14} It was shown that the generated friction is inversely proportional to tube separation, reaching values closer to the single-pore simulations (which should be considered as an infinitely-separated system) at the largest separation between the tubes, and estimation of an effective viscosity of the water molecules confined within the nanotubes confirmed this, along with consideration of the collective total velocity correlation functions (TCFs) for the water files within the arrays of CNTs.

A more general approach^{42,43} to the concept of the osmotic permeability based on linear response theory leads to the recognition that p_{f} is simply related to the slip-diffusion (or viscous slip) coefficient of a fluid undergoing pressure-driven transport in a nanopore structure. The results reported in this work also support the approximate relationship p_{f} ~ (N + 1)p_{d} suggested by Tajkhorshid et al ^{10} where N is the average number of water molecules within a single carbon nanotube.

30.MackerellA. D.Jr BashfordD.BellotM.DunbrackR. L.Jr EvanseckJ. D.FieldM. J.FischerS.GaoJ.GuoH.HaS.etal. J.PhysChem.102ashford, M. Bellot, R.L. Dunbrack. Jr., J.D Evanseck, M.J. Field, S. Fischer, J. Gao, H.Guo, S. Ha et al, J. Phys. Chem. B 102, 3586 (1998).

39.EnglishN. J.MacJ. M. D.ElroyJ.ChemPhys.119 11806 (2003).

40.EnglishN. J.MacJ. M. D.ElroyJ.ChemPhys.118 1589 (2003).

41.A.C. Metaxas and R.J. Meredith, ‘Industrial Microwave Heating’, Peter Peregrinus, London (1983).

42.S.SuhH.MacJ. M. D.ElroyMol.Phys58445 445 (1986).

43.MacJ. M. D.ElroyS.SuhH.MolPhys.60475 475 (1987).

44.Notethat.thepartition.coefficientK.relatesthe.equilibriumconcentrations.ofthe.waterbetween.thebulk.liquidvolume.themembrane.asvolumea.whole. C. N. T.poreslipidlayer.whereasKr.inEqn.6 is the equilibrium relation involving the CNT pore volume alone.

45.D.J. Evans and G.P. Morriss, ‘Statistical Mechanics of Non-equilibrium Liquids’,2 ed., Cambridge (2008).

46.Mc QuarrieD. A.Statistical‘.Mechanics’2 ed., University Science Books (2000).

Written By

Niall J. English, José-Antonio Garate and J. M. Don MacElroy

Submitted: April 6th, 2011Published: August 9th, 2011