Dynamic and weak electric double layers in ultrathin nanopores

Recent advances^1–3 in nanotechnology have led to the fabrication of nanopores in ultrathin membranes (e.g., single-layer graphene and single-layer MoS₂), opening up opportunities to explore a wide range of novel applications, such as water desalination,^4–9 power generation,^10,11 and disease diagnosis.^12,13 Of particular interest is the unique structural and dynamical properties of aqueous electrolytes in nanopores.¹⁴ Due to the unique chemistry and geometry of ultrathin nanopores, they have been shown to effectively reject salt ions, allowing water molecules to pass through the pores at high rates.^4,5,15 In addition, high ionic mobility combined with the strong charge selectivity of ultrathin nanopores has been shown¹⁰ to produce a remarkable nanopower density of up to 10⁶ W m⁻² when a salinity gradient exists across the membrane. The high ionic mobility in ultrathin nanopores is due to the weakened ionic adsorption to the thin surface of the pore.¹⁰ Ionic current blockade, which is used as a molecule identification signal in DNA sequencing and single biomolecule detection,^16–18 is also shown to possess a high resolution signal in ultrathin nanopores as the thickness of the ionic blockade region in the pore is comparable to the size of the DNA nucleotides and amino acids of proteins. The unique observations in ultrathin nanopores, from both experiments and computations, are mainly because of the distinctive structure and dynamics of molecules inside these ultrathin confinements as compared to those of thicker nanopores [e.g., carbon nanotubes (CNTs) and nanochannel slits]. In a thick nanopore, where the size of the solid surface is much larger than the size of ions, an electric double layer (EDL) is formed at the solid–liquid interface.¹⁹ Under the application of an electric field along the axis of the pore, EDLs in thick nanopores drag water molecules (solvent), resulting in a hydrodynamic flow known as electroosmosis, which has found numerous applications in micro- and nano-fluidic systems,²⁰ fuel cells,²¹ etc. In liquid dielectric-based field effect transistors (FETs), the presence of an EDL above the channel is shown to screen the charge of the target molecules (e.g., DNA molecules), which are outside the EDL away from the solid surface, hindering the sensitivity of the detection of target molecules.^22,23 The sensitivity increases with the increase in the thickness of the EDL characterized by the Debye length. Hwang et al.²³ showed that crumpling the graphene enhances the sensitivity of FETs by either increasing the Debye length in the deep concave surfaces or weakening the EDL above the sharp convex surfaces (which resemble thin edges of pores in ultrathin membranes). EDLs have also drawn a great deal of attention in the field of batteries and supercapacitors.^24,25 EDL capacitors are based on the adsorption of ions to the solid surface where a nanometer-long charge separation distance exists at the interface.²⁶ Given the fact that ultrathin nanoporous membranes are emerging as one of the most important materials for various applications, understanding the EDL or the ionic structure at their thin nanometer interface is of crucial importance.

Helmholtz²⁷ discovered that counter-ions in a solution are attracted to a charged surface (electrode) to neutralize the surface charge. In the Helmholtz model, the surface charge is counterbalanced only by the counter-ions, and the electric potential changes linearly from the surface to the counter-ions, leading to a constant capacitance. Helmholtz’s description of the EDL fails to account for the ions adsorbed onto the surface as well as the diffusing ions in the solution. In subsequent work, Gouy and Chapman²⁸ independently discovered that the capacitance changes with the application of an external potential, questioning the Helmholtz model. In the Gouy–Chapman model, the surface charge is neutralized by diffusing counter-ions in the solution, known as the diffuse layer, in which the variation of the ion concentration next to the surface obeys the Boltzmann distribution. The Gouy–Chapman model, however, fails in highly charged surfaces as it does not account for the finite size of ions. Stern²⁹ modified the Gouy–Chapman model by adding another layer, known as the Stern layer, where ions are ideally immobile and adsorbed onto a highly charged surface.

The interfacial properties in nanopores can be vastly different from those predicted by the classical models. As the thickness of the pore becomes comparable to the size of an ion in an ultrathin nanopore, the ionic structure is expected to differ from that of the classical EDL. When comparing the ionic structure in a finite-thickness nanopore to that of the classical EDL theory, one should note that the surface area is assumed to be infinitely large with an infinite number of particles in the classical theory. Here, using extensive molecular dynamics (MD) simulations, we study how the structure and dynamics of aqueous ions in single-layer graphene nanopores differ from those of thick CNTs. In graphene nanopores, the variation of the ion concentration and diffusion coefficient indicate that the ionic structure resembles mostly that of the bulk structure with no apparent Stern layer even for highly charged nanopore inner surfaces. In addition, ions and surface atoms are shown to be weakly correlated, resulting in a weak ionic adsorption onto the surface as compared to that in the thicker counterpart pores (e.g., CNTs). Due to the ultrathin nature of the single-layer graphene nanopore (0.34 nm thick), the residence time of a counter-ion, obtained from the correlations, is shown to be several orders of magnitude smaller than that of a thicker nanopore, resulting in the ultrafast dynamics of ions in ultrathin nanopores.

Each CNT simulation involves water molecules, ions (potassium and chloride ions), two graphene sheets (the entrance and exit walls), and a CNT between the sheets (Fig. 1). For the ultrathin graphene nanopore, there is only one graphene sheet with water molecules and ions [Fig. 1(a)]. (14,14) and (37,37) CNTs with thicknesses (h) of 3.34 and 10.34 nm are considered (the dimension in z, which is along the axis of the pore, varies from 6 to 16 nm depending on the thickness). In the case of graphene nanopores, a pore with a size similar to that of CNTs is drilled in the middle of the graphene sheet. As shown in Fig. 1(c), the graphene thickness is the physical size of a carbon atom defined by its Lennard-Jones (LJ) diameter (h = σ_c-c = 0.34 nm). To maintain the degree of porosity of the nanopores, the dimensions in x and y (which lie in the plane of graphene sheets) vary from 6.25 to 12 nm depending on the diameter (see the Methods section for more details). In the equilibrated MD simulations, the number of potassium ions (K⁺, which are the counter-ions as the carbon atoms of the pore inner surface are assigned negative charges with charge densities ranging from −0.0006 to −0.6 C m⁻²) is calculated inside an ultrathin single-layer nanopore (h = 0.34 nm) and a thick CNT (h = 10.34 nm) with a diameter of 1.9 nm. Assigning charges to the carbon atoms of a graphene nanopore may seem physically questionable; however, it has been shown in experiments¹¹ that ultrathin graphene nanopores exhibit charge selectivity upon changing the pH of the solution. This is due to the existence of charged functional groups at the edge of the pore. In our simulations, the effect of these charges is mimicked by simply assigning charges to the carbon atoms at the edge of the pores. The number of K⁺ ions, which is an indication of the presence of an ionic layer or EDL inside the nanopores, is shown as a function of the simulation time in Fig. 2. In the ultrathin nanopore, the presence of K⁺ ions has a dynamic nature where an ionic layer exists intermittently in the nanopore with a time scale of the order of a few picoseconds. As the thickness of the nanopore increases slightly to h = 3.34 nm, the ionic structure inside the nanopore begins to assume a continuous and rigid EDL, which is stable over the entire simulation time (>nanoseconds).

The nature of the ionic structure inside the nanopores changes abruptly as we approach the ultrathin limit of a single-layer graphene. To shed more light on the structure and dynamics of K⁺ counter-ions inside the nanopores, the ion concentration and the axial (z-direction) diffusion coefficient are plotted as a function of the radial distance from the center of the pore in Fig. 3. The diffusion coefficients in the z-direction are calculated locally within the physical thickness of the pore using the mean-squared displacement (MSD) $Dz=zit−zi022t$⁠, where z_i is the position of ion i at time t. Locally, the MSD of an ion entering the pore region is calculated while the ion remains inside the pore region (the bin of interest). An averaged MSD is then obtained for all instances of ions entering the pore region during the simulation time. In the case of large-thickness CNTs (h = 10.34 nm), an EDL is formed for both low and high surface charge densities. The structure of the EDL in the low surface charge density CNT is consistent with the structures obtained in previous studies,^19,30 where a counter-ion concentration peak is formed next to the wall and the co-ion concentration is significantly reduced near the wall (with an ionic charge density of 0.207e nm⁻³). The counter-ion diffusion coefficient along the axis of the pore is almost uniform with no features of a Stern layer. However, for the high surface charge density, we observe a large counter-ion concentration peak (with a high ionic charge density of 2.032e nm⁻³) and an almost vanishing ion diffusion near the charged wall. In addition, a strong layering of both co- and counter-ions is formed within just 0.2 nm of the wall surface. Comparing the ionic structure and dynamics obtained from MD simulations in thick CNTs to the ones predicted by the classical EDL models, the EDL in the low surface charge density CNTs resembles the Gouy–Chapman model, where a diffuse layer with a higher concentration of counter-ions is present [Fig. 4(a)]. The ion concentration variation, however, does not follow the Boltzmann distribution as assumed in the Gouy–Chapman model. The properties of the EDL in large surface charge density CNTs are similar to those of the Gouy–Chapman–Stern model [Fig. 4(b)], where an additional immobile Stern layer is formed next to the surface. In the supplementary material, an animation is provided where the immobility of ions is clearly displayed. For the nanopore in a single-layer graphene with a low inner pore surface charge density [Fig. 3(c)], the ion concentrations (with an almost neutral ionic charge density of −0.008e nm⁻³) are nearly homogeneous throughout the pore with a lower concentration near the surface. The region next to the surface is highly concentrated with water molecules. The counter-ion diffusion coefficient is nearly uniform. The structure and dynamics of ions inside the low surface charge density graphene pore do not resemble those of the existing EDL models as neither a diffuse layer nor a Stern layer is present [see Fig. 4(c) for the schematic of ions inside low surface charge density nanopores]. Rather, despite the nanopore having a net charge, the nanopore acts like a point of zero charge.³¹ In the high pore surface charge density graphene nanopore [Fig. 3(d)], near the wall (which is highly concentrated with water molecules), a weak and dynamic diffuse layer [Fig. 4(d)] is formed where the concentration of counter-ions is lower or comparable to the bulk concentration (with a small ionic charge density of 0.079e nm⁻³). This diffuse layer differs from that of the Gouy–Chapman (GC) model as the counter-ion concentration in the GC model is expected to increase (higher than the bulk value) near the wall. As the diffusion coefficient is nonzero near the wall, there is no indication of a Stern layer even for such a high pore surface charge density ultrathin nanopore.

To further investigate the properties of the EDLs and ionic layers in ultrathin graphene nanopores and CNTs, the charge screening factor of ions is studied in Fig. 5. The screening factor of the surface charge by the ions is defined by $SF=∫0rq(r)drσ$⁠, where r is the radial distance away from the pore surface and q(r) is the net ionic charge density. As shown in Fig. 5, for the low surface charge density CNT, the surface charge is screened within 1 nm of the CNT surface. For the high surface charge density CNT, the surface charge is over-screened within a short distance from the surface (∼0.3 nm) due to the presence of the Stern layer. This over-screening and charge inversion³² cannot be explained by the classical EDL theories. For the graphene nanopores, in both low and high surface charge densities, the surface charge is never completely screened by the ions inside the pores (only a fraction of the surface charge is screened). The surface charge must therefore be compensated at large distances outside the pore. The weak ionic screening in graphene nanopores is due to the lack of a high enough counter-ion concentration near the pore walls. To study the residence times of the counter-ions near the wall (the degree to which the counter-ions are adsorbed onto the wall surface), the time correlation¹⁹ between the counter-ions and the wall atoms is calculated using $Cijt=θij0θijt$⁠, where i and j are a pair of particles (here, the K⁺ ion and carbon atom are the particles). θ_ij is 1 if i and j are the nearest neighbors if they are within the first peak in the concentration plots and 0 otherwise. C_ij is the probability of a pair of particles being the nearest neighbors from time 0 to t; therefore, the residence time of the K⁺ ion near the wall can be obtained from $τ=∫0∞Cijtdt$⁠. In Fig. 5, the K⁺-wall time correlation is plotted as a function of time in the ultrathin graphene nanopore and 10.34 nm-thick CNT for different surface charge densities. In the graphene nanopore, the correlations die out rapidly regardless of the magnitude of the surface charge density. The calculated K⁺ ion residence times are about 4.85 ps for the various surface charge densities considered. The short residence time suggests that a Stern layer is not present in ultrathin graphene nanopores for the practical range of surface charge densities considered in this work. For the thick CNTs, the correlations die out quickly for low surface charge densities with K⁺ ion residence times of about 15 ps. However, for large surface charge densities, K⁺ ions and wall carbon atoms are highly correlated with much larger residence times (beyond the time scales easily accessible with MD simulations). To estimate the residence times for CNTs with large surface charge densities, the tail of C_ij(t) is extrapolated by fitting a curve to the data from MD simulations. The residence times based on the extrapolated C_ij(t) are found to be of the order of seconds. The strong K⁺-wall correlation is a signature of ion immobility or the presence of a Stern layer near the wall. Therefore, as the thickness increases from the ultrathin limit of a single-layer graphene to thicker CNTs, the residence times (or the degree of ionic adsorption onto the surface) increase abruptly by about 12 orders of magnitude. In the ultrathin limit of a graphene nanopore, no Stern layer is observed regardless of the strength of the surface charge density, which justifies the fast dynamics of ions inside these nanopores that has been observed both computationally and experimentally.¹⁰ 

The structural and dynamical properties of aqueous ions in ultrathin nanopores are vastly different from those predicted by the classical theories. MD simulations indicate that the EDL at the nanopore interface is weak and exhibits a dynamic behavior as the thickness of the nanopore reduces to that of a single-layer graphene nanopore. Unlike the EDL in thick CNTs, where diffuse and Stern layers are observed, in ultrathin nanopores, no Stern layer is present, and a weak diffuse layer is present only for highly charged nanopores. The lack of a Stern layer formation is shown to be due to the very weak correlation between the counter-ions in the interfacial layer and the wall atoms. The weak correlation (or weak ionic adsorption) results in fast dynamics of ions with short residence times in ultrathin graphene nanopores. Analyzing the correlation functions reveals that a huge reduction (from seconds to picoseconds) is observed in residence times as the thickness decreases to the limit of a single-layer graphene nanopore. The critical thinness is reached only for single-layer graphene nanopores (and other ultrathin 2D materials such as MoS₂) for which ultra-weak correlation (adsorption) and ultra-fast residence times of counter-ions are observed.

MD simulations were performed using the Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) package.³³ The simulation boxes contain 20 000–100 000 atoms depending on the thickness and diameter of pores with ions (potassium and chloride) having a molarity of ∼1M. The Extended Simple Point Charge (SPC/E) water model was used, and the SHAKE algorithm was employed to maintain the rigidity of the water molecule. For nonbonded interactions involving ion–ion and ions with water and carbon atoms, the mixing rule was employed to obtain the LJ parameters.³⁴ The carbon–water interactions were modeled by the force-field parameters given by Wu and Aluru.³⁵ Since the carbon atoms are fixed in space, the interactions between carbon atoms were turned off. The effect of flexibility is neglected for simplicity and the fact that fixing carbon atoms does not lead to a significant change in the structure and dynamics of fluids (e.g., flow rate errors of about 10% are observed when using fixed carbon atoms in CNTs).^14,36 The LJ cutoff distance was set to 12 Å. The long-range electrostatic interactions were calculated by the Particle–Particle–Particle–Mesh (PPPM) method.³⁷ Periodic boundary conditions were applied in all three directions. For each simulation, first, the energy of the system was minimized for 10 000 steps. Next, the solution was equilibrated using the NPT ensemble for 1 ns at a pressure of 1 atm and a temperature of 300 K with a time step of 1 fs where the simulation size fluctuations are only allowed along the axis of the pore. With the graphene atoms held fixed in space, the NPT simulations allow the water to reach its equilibrium density (1 g/cm³). Then, an additional NVT simulation was performed for 2 ns to further equilibrate the system. Temperature was maintained at 300 K by using the Nosè–Hoover thermostat with a time constant of 0.1 ps.^38,39 Finally, the production equilibrium simulations were carried out in the NVT ensemble for ∼20 ns.

See the supplementary material for the animations of the trajectories of ions inside the pores considered in this work.

The research on CNTs was supported as part of the Center for Enhanced Nanofluidic Transport (CENT), an Energy Frontier Research Center funded by the U.S. Department of Energy, Office of Science, Basic Energy Sciences under Award No. DE-SC0019112. The other aspects of the work presented here were supported by the National Science Foundation under Grant Nos. 1545907, 1708852, 1720633, and 1921578. The simulations were performed using the Extreme Science and Engineering Discovery Environment (XSEDE) [supported by National Science Foundation (NSF) Grant No. OCI1053575] and Blue Waters (supported by NSF Award Nos. OCI-0725070 and ACI-1238993 and the state of Illinois).

The data that support the findings of this study are available within the article and its supplementary material and from the corresponding author upon request.