Optimal transport and colossal ionic mechano-conductance in graphene crown ethers

TL;DR

This study demonstrates that applying small mechanical strains to graphene crown ether pores dramatically enhances ion conductance and selectivity, providing a tunable platform for investigating optimal ion transport mechanisms.

Contribution

The paper introduces a novel graphene-based platform with mechanically tunable ion transport properties, revealing colossal conductance changes and insights into biological ion channel function.

Findings

Small pore strains induce 100% increase in conductance.

Transport operates mainly in a diffusive, barrierless regime.

Mechanical modulation reveals local electrostatic conditions.

Abstract

Biological ion channels balance electrostatic and dehydration effects to yield large ion selectivities alongside high transport rates. These macromolecular systems are often interrogated through point mutations of their pore domain, limiting the scope of mechanistic studies. In contrast, we demonstrate that graphene crown ether pores afford a simple platform to directly investigate optimal ion transport conditions, i.e., maximum current densities and selectivity. Crown ethers are known for selective ion adsorption. When embedded in graphene, however, transport rates lie below the drift-diffusion limit. We show that small pore strains -- 1 % -- give rise to a colossal -- 100 % -- change in conductance. This process is electromechanically tunable, with optimal transport in a primarily diffusive regime, tending toward barrierless transport, as opposed to a knock-on mechanism. Measurements of mechanical current modulation will yield direct information on the local electrostatic conditions of the pore. These observations suggest a novel setup for nanofluidic devices while giving insight into the physical foundation of evolutionarily--optimized ion transport in biological pores.

Equations

$$\Delta F_{\mathrm{K}}=\Delta E_{\mathrm{deh}}+E_{\mathrm{KO}}+E_{\mathrm{KC}}+E_{\mathrm{KK}},$$

$$\chi =\frac{d\Delta F}{da}=\sum _{⟨\mu \nu ⟩}\frac{\partial \Delta F}{\partial {\rho }_{\mu \nu }}\frac{d{\rho }_{\mu \nu }}{da}\approx q_{\mathrm{K}}{E}_{\rho },$$

$$I\approx A{e}^{-\Delta F^{\star }/k_{B}T},$$

$$\frac{\Delta I}{I}\approx e^{(-\Delta F^{\star }(s)+\Delta F^{\star }(s=0))/k_{B}T}-1\approx e^{\alpha as\chi /k_{B}T}-1,$$

$$I_{\mathrm{K}}$$

Full text

††volume: 00††issue: 0

Introduction

The search for universal transport mechanisms is an underlying theme in ion channel research (1). Nonetheless, even foundational questions – such as the primary mechanism for selectivity of K+ over Na+ in the potassium ion channel family – remain contentious: Is the smaller dehydration energy of the K+ ion responsible, or it is the “snuggle–fit” of K+ in the selectivity filter (2, 3, 4)? Solving these puzzles requires a detailed means to control and characterize ion transport. Biological ion channels are notoriously difficult to manipulate, requiring mutagenesis and near–native conditions for patch–clamp studies, thus limiting experimental characterization and tunability. These large macromolecules are also theoretically complex, requiring hefty resources to conduct statistically meaningful simulations. In contrast, biomimetic synthetic pores – or, as addressed here, pores that retain core aspects of some biological channels – may provide a simple platform to explore ion transport mechanisms under a broad range of conditions (5).

In this work, we demonstrate that graphene crown ether pores (6) have competing electrostatic and dehydration contributions to transport, reminiscent of biological ion channel mechanism and selectivity. The discovery of crown ethers a half–century ago triggered a cascade of research that evolved into the fields of host–guest and supramolecular chemistry (7). The isolation of graphene (8) has likewise inspired diverse fields of study, including those that leverage atomically thin pores for biosensing and nanoscale separation (9). Graphene crown ether pores lie at the crossroads, affording a means to address transport and selectivity in a synthetic analog of biological ion channels. Using a graphene–embedded 18-crown-6 ether pore as a model system, we examine the mechanism for ion translocation under variable strain and applied voltage, see Fig. 1. Our results indicate that a few percent change in pore size can induce a several fold increase in ionic current, underscoring the sub–nanometer fine tuning between electrostatics and dehydration seen in biological contexts. We similarly find that small changes in pore radius and charge can optimize transport, and even alter the underlying mechanism of ion translocation. These results collectively suggest that graphene crown ether pores are an effective mechanistic probe for understanding transport in both biological and artificial ion channels alike.

Results

Partial charge distribution: Unlike an isolated crown ether, an 18-crown-6 pore in graphene is believed to assume a planar conformation (6). However, while critical to molecular dynamics (MD) simulations, the charge distribution and response of this pore is unknown. Prior simulations on this and related systems assume oxygen partial charges that span between $q_{\mathrm{O}}=-0.21\,e$ and $q_{\mathrm{O}}=-0.74\,e$ (10, 11, 12, 13, 14, 15), over which transport properties will vary dramatically. To provide a more comprehensive approach, we reexamine this assignment using density functional theory (DFT), see the Methods and supplementary materials (SM) for details. Determinations from electrostatic potential fitting without an ion present – consistent with additive force fields and the electrostatic environment predicted by DFT – yield a uniform value of $q_{\mathrm{O}}=-0.24\,e$, decreasing to $-0.23\,e$ when a K+ is present. Chemically distinct, but oxygen containing, fragments in OPLS and CHARMM force fields give substantially higher partial charges, $\geq-0.4\,e$, of which we take $q_{\mathrm{O}}=-0.54\,e$ as a representative example. An upper bound of $q_{\mathrm{O}}=-1.0\,e$ is given by Bader analysis, which tends toward larger values of partial charge and represents the maximum expected from any electronic structure calculation. Physically reasonable values of $q_{\mathrm{O}}$ should span the low end of the range $-0.2\,e$ to $-0.6\,e$, as evidenced by electrostatic fits of $q_{\mathrm{O}}=-0.34\,e$ for the 18-crown-6 alkyl ether and $q_{\mathrm{O}}=-0.54\,e$ for a highly–dipolar glycine dipeptide backbone carbonyl. The graphene crown ether should be lower than both. We perform MD simulations using three distinct assignments $q_{\mathrm{O}}=(-0.24,-0.54,-1.0)\,e$ to capture the full range of potential behavior, showing both universal features of transport and aspects particular to different assignments. We test the effect of strain on ion transport using different applied biases and deformations, well within the elastic limit of graphene (16, 17).

Ion transport mechanism: Since this pore is tiny and negatively charged, only K+ contributes to the current in a KCl solution (set at a concentration of 1 mol/L in our case). The mechanism of transport depends on the ion’s binding strength in the pore, largely regulated by $q_{\mathrm{O}}$ (as an indicator of the oxygen–carbon charge separation) and atomic structure. The mechanism is normal drift–diffusion for $q_{\mathrm{O}}=-0.24\,e$, shifting to a knock–on type for $-0.54\,e$ and $-1.0\,e$ (see the illustration in Fig. 1). For drift-diffusion, the residence time (on the order of 0.1 ns at 0.25 V) is much shorter than the delay to the next ion crossing event (on the order of 1 ns at 0.25 V). Conversely, for the two-step knock-on mechanism seen at higher charge, the residence time (on the order of 1 ns) is much larger than the delay time (less than 10 ps for most events; see the SM). However, even for $q_{\mathrm{O}}=-0.54\,e$, drift–diffusion starts to dominate at high strain and/or voltage, both of which weaken ion binding to the pore.

Ionic mechano-conductance: Figure 2a depicts a potential experimental setup, in which a graphene sheet is selectively strained by deforming an underlying substrate. Our simulations mimic these conditions by increasing the cross–sectional area of the membrane to induce a strain. The relative change in current as a function of strain, for different values of $q_{\mathrm{O}}$ and applied voltages, is shown in Fig. 2b. At low voltages (0.25 V and 0.5 V), the ionic current increases rapidly with strain, orders of magnitude larger than expected from the increase in pore area alone. Instead, this colossal increase is indicative of a change in the energetics and dynamics of ion transport. In the large voltage regime (1 V), the current remains constant with increasing strain when $q_{\mathrm{O}}=-0.54\,e$. This is not altogether surprising, as a large voltage drop across the membrane suppresses the influence of the free energy profile on translocating ions. The dependence of the current on $q_{\mathrm{O}}$ and $V$ is in Fig. 2c. In this case, the current increases with $q_{\mathrm{O}}$ between $q_{\mathrm{O}}=-0.24\,e$ to $q_{\mathrm{O}}=-0.54\,e$, while giving a turnover to decreasing current at higher charge.

The energy landscape: Our results indicate that the energetics of pore permeation lie in a regime where small electromechanical variations modify current magnitudes and give rise to qualitative changes in behavior. The free energy of a potassium ion in the vicinity of the pore, $\Delta F_{\mathrm{K}}$, is approximately the sum of a dehydration energy barrier $\Delta E_{\mathrm{deh}}$ and three electrostatic terms ($E_{\mathrm{KO}}$, $E_{\mathrm{KC}}$, and $E_{\mathrm{KK}}$) that capture the interaction of K+ with nearby charged atoms:

[TABLE]

where $E_{\mathrm{deh}}=\eta\sum_{i}f_{i}E_{i}$ is the barrier due to dehydration, $f_{i}$ ($E_{i}$) is the fractional dehydration (energy) of the $i^{\mathrm{th}}$ hydration layer, and the factor $\eta\approx 1/2$ accounts for the nonlinear effects  (18). Each electrostatic component is $E_{\mu\nu}=q_{\mu}n_{\nu}q_{\nu}/4\pi\epsilon_{0}\epsilon(z_{\mu\nu})r_{\mu\nu}$, where $\epsilon(z_{\mu\nu})$ is a position–dependent relative permittivity (19), $q_{\mu(\nu)}$ is the charge of ion species $\mu(\nu)$, $n_{\nu}$ is the number of proximate species $\nu$ (6 for oxygen, 12 for carbon, and 1 or 0 for a nearby K+), $z_{\mu\nu}$ ($r_{\mu\nu}=\sqrt{z_{\mu\nu}^{2}+\rho_{\mu\nu}^{2}}$) is the axial (total) distance between species $\mu$ and $\nu$, and $\rho_{\mu\nu}$ is the respective radial distance. $E_{\mathrm{KK}}$ is the energetic coupling of an incoming K+ to a trapped K+. The van der Waals interactions, the coupling to other solvated ions, and the entropy change should have only marginal contributions when $\Delta F_{\mathrm{K}}$ varies with strain (e.g., the entropy change into the pore is important, but this varies little with moderate strain). In order to quantify the variable electrostatic environment reflected in $\epsilon(z)$, as well as nonlinear effects in $\Delta E_{\mathrm{deh}}$ (stronger polarization of bound water molecules as an ion becomes dehydrated (18)), we employ all–atom MD simulations to find the free energy profile of K+ around the pore (Fig. 3 and the Methods). Equation 1 explains the key features in the free energy profiles for the range of $q_{\mathrm{O}}$ under consideration. Together with the MD data, this expression dissects the contributions to the free energy. We will use this equation below to characterize the response of ion transport to strain.

Figure 3 shows that when $q_{\mathrm{O}}=-0.24\,e$, a free energy barrier appears outside the pore due to the progressive dehydration of K+. The barrier decreases with increasing proximity to the pore opening as the flanking oxygen atoms begin to compensate for water molecules displaced from the ion’s first hydration layer. This results in a potential well on both sides of the membrane. Further dehydration occurs as the ion enters the pore, giving rise to an overall barrier at $z=0$ since the electrostatic attraction between the pore oxygen atoms and the K+ ($E_{\mathrm{KO}}$) is insufficient to compensate for the dehydration energy penalty. This central barrier – and only shallow adjacent wells – ensure that an approaching K+ will almost always find an empty pore, and hence the $E_{\mathrm{KK}}$ term will not contribute in this case. Equation 1 confirms this assignment of the peaks and wells.

This situation is inverted for $q_{\mathrm{O}}=-0.54\,e$, see Fig. 3, where a deep potential well forms at the center of the pore due to strong electrostatic attraction with the K+. The central pore acts as a binding site in this case, ensuring that it is almost always occupied. When K+ approaches the pore, the presence of an interstitial K+ repels it, yielding a potential barrier on either side of the pore. A free energy calculation with single K+ and single Cl- confirms this (see the SM). We can estimate $\epsilon_{r}(z)$ using the free energy profile from Fig. 3 and the fractional dehydration from MD. For the physical range of pore charge, $\epsilon_{r}(z)$ is around 4 to 6 within 0.2 nm of the pore, which is consistent with experimental values for the dielectric constant within 1 nm of an interface (20). A similar value for dielectric constant is given in another recent computational study (19).

Barrierless transport: While the transport mechanism is different for unstrained pores, the colossal mechano–conductance change is present at low voltage regardless of the partial charge. In both partial charge scenarios, the free energy profile becomes smoother with increasing strain, tending toward “barrierless transport” (21, 19): Pulling charged oxygens away from the pore center reduces their coupling to translocating ions. The colossal change in conductance indicates that there are optimal structural positions for the oxygens: Picometer changes in their position do not change the dehydration contribution to the central barrier (see the SM) but it vastly changes the transport rate, maintaining exclusion of other ions but enhancing transport rates by several fold. Biological ion channels can further make use of the partial charge, engineering not only structural characteristics, but also the electronic environment. As we discussed earlier, the larger partial charges we consider – in line with biological channels (22) – show a shift from knock–on to drift–diffusion mechanisms as transport starts to become barrierless. In this manner, the electromechanical environment may exert a large effect on the current.

Mechanical susceptibility: To understand the mechanics underlying this process, we define the susceptibility, $\chi$, of the free energy to a small change in pore size,

[TABLE]

where $a$ is the effective pore radius, $\langle\mu\nu\rangle$ indicates the pairs KO and KC (we ignore the role of KK interactions, although this is important for knock–on aspects of transport), and $d\rho_{\mu\nu}/da\approx 1$. The latter approximation is reasonable since, as the pore size increases, the oxygen and carbons move outward by the same amount. For simplicity, we took this for an ion at the origin ($z=0,\rho=0$) where the susceptibility is just the local radial electric field, $E_{\rho}$, multiplied by the charge of the translocating ion. Due to the lack of electrostatic screening in the pore (i.e., $\epsilon(z=0)$ is small) and the proximity of the oxygen atoms, this field is very large, resulting in a large $\chi$: Equation (1) gives $\chi\approx 100\,k_{B}T$/nm and $\chi\approx 140\,k_{B}T$/nm for the $q_{\mathrm{O}}=-0.24\,e$ and $q_{\mathrm{O}}=-0.54\,e$, respectively (using $\epsilon=3.9$ and $6.4$, which we extract via Eq. (1) and the all-atom MD data). These two susceptibilities are closer than their difference in charge would indicate due to the smaller permittivity for $q_{\mathrm{O}}=-0.24\,e$ (the charge $q_{\mathrm{O}}=-0.54\,e$ more strongly attracts counter charges, specifically K+ and hydrogens on water molecules, that help screen the interaction.)

The large susceptibility means that even a one picometer change in the pore size gives a 0.1 $k_{B}T$ to 0.14 $k_{B}T$ change in the free energy. A 1 % strain gives about 7 pm change in radius. Using this in Eq. (1), free energy changes by about 0.7 $k_{B}T$ for $q_{\mathrm{O}}=-0.24\,e$ and about 1.0 $k_{B}T$ for $q_{\mathrm{O}}=-0.54\,e$, in line with the all-atom results in Fig. 3. We note that there are clearly other effects occurring: As we will discuss below, the $q_{\mathrm{O}}=-0.24\,e$ case is already in a limit where dehydration and electrostatic effects are comparable, and thus the change in the free energy profile versus $z$ reflects both (the barrier reduction at $z\approx 2$ nm is due to a changing dehydration contribution, whereas in between those two outer barriers, electrostatic effects dominate and $\Delta F$ increases overall). For $q_{\mathrm{O}}=-0.54\,e$, the potential well shallows, causing a transition from a knock–on to a diffusive mechanism and removing the ion–ion repulsion that causes the external barrier. Nevertheless, in both cases, the electrostatic contribution to the susceptibility captures the increased free energy for an ion in the pore. While a general susceptibility will need to include the totality of effects, electrostatics and dehydration will typically dominate and are key to understanding amplification.

The free energy profile and its variation (according to $\chi$) determines the current. In particular, the largest free energy hurdle, $\Delta F^{\star}$, in the profiles of Fig. 3 are the main factors influencing the conductance. For $q_{\mathrm{O}}=-0.24\,e$, this will be the jump from about $z=$ 0.1 nm to 0.2 nm and, for $q_{\mathrm{O}}=-0.54\,e$, this will be the jump out of the well ($z=0$) to the maxima (at about $z=0.2$ nm to $z=0.4$ nm depending on strain). The current will take the form (see Methods)

[TABLE]

where $A$ is a constant that only weakly depends on pore parameters (such as size), but we will take it as fixed. The total amplification of the current then becomes

[TABLE]

where $-\Delta F^{\star}(s)+\Delta F^{\star}(s=0)\approx d\Delta F/da\cdot\delta a\approx\alpha\,a\,\,s\,\chi$. The parameter $\alpha=1/a\cdot da/ds$ quantifies the relative response of the pore radius to the application of a strain. Both MD and DFT predict that $\alpha\approx 2$, meaning that if the graphene has 1 % strain, the atoms at the pore rim move by about 2 %. This is due to a distortion of the rim structure since the pore is an effective impurity and can relieve strain by relaxing angular coordinates. The straightforward relation in Eq. 4 indicates that the gain in current is the exponential of the work performed by the applied strain on a K+ within the pore. For the values in our case, Eq. (4) gives an amplification of 100 % to 200 % for a 1 % strain, in line with the MD results in Fig. 2b.

Optimum transport: In all parameter regimes, the amplified current is due to an overall flattening of the free energy profile. Ion transport thus tends toward barrierless transport for small values of strain and displays a colossal increase in conductance along the way. However, for $q_{\mathrm{O}}=-0.24\,e$, the dominant contribution to $\Delta F^{\star}$ shifts from a pathway connecting a minimum (at $z=\pm\,0.1$ nm) and the external solution to a pathway that hops over the central barrier (the well minima already lie at about $\Delta F=0$). After this point, $\Delta F^{\star}$ will continually increase due to decreasing compensation by the oxygen partial charge. There will be a corresponding decrease in current, see Fig. 4, thus indicating an optimum near this value of strain – a regime that may be experimentally accessible.

In separate calculations with an equal mixture of NaCl and KCl (with Cl- at one molar), we find that only a single Na+ traversed the pore (compare to 47 K+) during 1 $\mu$s of simulation for the unstrained pore at a 0.25 V bias. At higher strains, only K+ translocate through the pore within the 1 $\mu$s simulation time. The lack of Na+ transport is due to the larger hydration energy of Na+ compared to K+, giving it a higher penalty to go through the pore. The presence of an overall barrier should result in Na+ crossings following Poisson statistics. Due to rarity of these crossings, we consider two values for the Na+ current, 0.5 pA and 0.16 pA . The former would give a 5 % probability to detect no Na+ crossings during 1 $\mu$s and thus is a high confidence upper bound to the Na+ current. The latter gives equal probabilities for 0 and 1 Na+ crossing and thus represents an estimate of the Na+ current assuming we are near the threshold to start seeing Na+ crossings, for which the single crossing at 0 % strain suggests is the case. Using the upper bound on the Na+ current, the corresponding high-confidence lower bound on the selectivity ranges from about 15 at 0 % strain to 52 at 3 % strain (note that the K+ currents in the equal mixtures are 7.5 pA at 0 % strain and 26.2 pA at 3 % strain, which are approximately half the value in the 1 mol/L KCl case shown in Fig. 4 due to the mixture having half the K+ concentration). The estimate of selectivity, though, ranges from 47 at 0 % strain to 164 at 3 % strain. Very long simulations (much greater than 10 $\mu$s) will be necessary to tightly bound the selectivity, but it will be large nonetheless.

For $q_{\mathrm{O}}=-0.54\,e$, the free energy at $z=0$ will not turn from a well to barrier until about $s=6$ % (i.e., about $\delta a=6\,k_{B}T/\chi$), with the optimum appearing at potentially larger strains and outside of measurable regimes. Up to this point, the current will increase with strain as the free energy profile becomes truly barrierless. In either case, the dehydration contribution itself (at $z=0$) will not substantially change until the pore becomes sizable, as indicated by the continuing decrease of the current in Fig. 4. The selectivity filters of biological channels have partial charges near $q_{\mathrm{O}}=-0.54\,e$ by virtue of the dipolar peptide backbone (4, 23, 24, 22), lying within this more strongly coupled regime and affording them extensive mechanistic leeway. These systems generically exploit the spatial distribution of charges or functional groups via the protein structure to achieve optimality. Strain in the graphene crown ether pore modulates a single collective variable, acting as a model for this optimization process.

Discussion

The investigation of 18-crown-6 and other crown ether pores, as well as pores generally, with respect to strain (the experimental extraction of $\chi$ and determination of the local interaction environment), voltage, and ion concentration will shed further light into the behavior of selectivity and optimal transport. These characteristics underlie the mechanism of physiological ion channels and are of substantial interest in biology and medicine (1, 25). This is paralleled by broad applicability to technologies for the sensing and sequencing of biomolecules (26, 27, 28); molecular sieving and desalination (29, 30); as well as batteries, fuel cells, and other energy harvesting devices (31, 32). In all of these applications, the atomic thickness and physical strength of graphene provide a unique advantage (28, 27) facilitating the detailed examination of molecular processes such ion dehydration (18, 33). Mechanical modulation of ion transport will thus lend itself to filtration and other technologies.

Our results indicate the ideal pore characteristics for optimal transport: Minute changes in pore size or charge give rise to a dramatic change in ionic conductance as transport becomes barrierless, but further change will decrease the conductance. These systems also comprise a simple and transparent platform to understand the evolution of biological ion channels. Working by analogy, it is possible to identify structural features and charge distributions that maximize selectivity while maintaining high permeation rates. These considerations, in turn, constrain which conformational transitions may be leveraged to regulate transport (34). Layered crown ether graphene pore heterostructures – mimicking aspects of the selectivity filter in the potassium channel from Streptomyces lividans (KcsA) – are likewise promising for workable models of biomolecular channel dynamics. These systems should be quite sensitive to ionic conditions, potentially gating or doping graphene to access broad regimes of behavior. The susceptibility to transverse strain is also generally interesting for artificial porous membranes as it gives otherwise inaccessible information regarding the local electrostatic and hydration conditions. This platform will thus open new and general opportunities for the study, quantification, and tuning of ion transport at the nanoscale.

Methods

We determine point charges using both localized–basis and extended plane–wave DFT schemes. The localized basis calculations (B3LYP/6-31G*) employ a truncated model of the graphene pore alongside the Gaussian09 code (35), affording point charges through CHELPG electrostatic potential fitting (36). In contrast, for plane–wave DFT, we employ VASP(37) and a dispersion–corrected scheme (PBE-D2/PAW) for a 1.71 nm $\times$ 1.73 nm graphene supercell containing an 18-crown-6 pore, with periodic images separated by a 2.0 nm vacuum layer ($5\times 5\times 1$ Monkhorst–Pack mesh, 500 eV cutoff) (38). In this case, point charges were determined using Bader analysis (39). A complementary set of plane wave calculations employs a similar pore within a 2.96 nm $\times$ 2.99 nm cell, described using CP2K (40) and a mixed Gaussian plane–wave basis scheme (PBE-D3/GPW/DZVP/GTH pseudopotentials, 5442 eV real–space grid cutoff) with $\Gamma$-point sampling and RESP charge fitting (41). Point charges are derived from optimized geometries, defined where SCF energy convergence is below $1.0\times 10^{-4}$ eV and forces below 0.1 eV nm*-1*. CHELPG fitting yields $q_{\mathrm{O}}=-0.24\,e$ without a K+ present in the pore and $-0.23\,e$ when a K+ is present. The complementary plane-wave calculations confirm these assignments. The Bader analysis yields $q_{\mathrm{O}}=-1\,e$, an upper bound for local partial charges due to its strong bias toward chemically intuitive, but not electrostatically representative, distributions. Both MD (below) and DFT give an $\alpha$ of about 2 (DFT yields a slightly higher value), indicating that the pore rim has an enhanced response to the strain. Reference fragments of $q_{\mathrm{O}}=-0.54\,e$ include calculated CHELPG carbonyl charges for the Ac-Gly-Gly-CH3 peptide (B3LYP/6-31G*), as well as standard OPLS ($q_{\mathrm{O}}=-0.50\,e$) and CHARMM ($q_{\mathrm{O}}=-0.51\,e$) backbone carbonyl parameters. Graphene should have oxygen partial charges lower than any of these fragments’ values – as well as the isolated crown ether’s oxygen partial charges – due to the delocalized nature of the electrons in the carbon sheet (see the SM for details and extended methods).

We simulate ion transport through the pore using all–atom molecular dynamics via the NAMD2 code (42) with an integration time step of 1 fs. We take a simulation cell aspect ratio of 1.2 – the golden aspect ratio – where the simulation captures both access and pore resistances (43, 44). While access resistance should be less than a 10 % correction for some parameter regimes, it must be taken into account when understanding how the drift–diffusion limit is approached. Its exclusion by choosing inappropriate (either too small or a poor aspect ratio) simulation cells hinders the comparison with the drift–diffusion limit (also see the SM).

We employ the adaptive biasing force (ABF) method (45) to calculate the equilibrium free energy $\Delta F$ profile for a K+ going through the pore. We compute $\Delta F$ versus $z$ in a cylindrical region of radius 0.1 nm and length 3 nm centered around the pore. As a K+ approaches the pore, it experiences either a free energy barrier or a potential well. Assuming a small applied voltage $V$, the current over a barrier $\Delta F>0$ is (46)

[TABLE]

and the current across a potential well $\Delta F<0$ is

[TABLE]

where $k_{\text{in}}$ and $k_{\text{out}}$ are the rate constant for incoming and outgoing ions and $c$ is the ion concentration in bulk solvent. Both of these expressions give Eq. 3 for large $\Delta F$. (See the SM for a additional discussion and results regarding the influence of the free energy barrier, including a multiple barrier rate model).

Acknowledgments

S. Sahu, J. Elenewski, and C. Rohmann acknowledge support under the Cooperative Research Agreement between the University of Maryland and the National Institute of Standards and Technology Center for Nanoscale Science and Technology, Award 70NANB14H209, through the University of Maryland.

Author contributions: S.S. analyzed the ion transport behavior and performed MD simulations. J.E. and C.R. performed DFT calculations. S.S., J.E., and M.Z. formulated and derived mathematical expressions. All authors developed the ideas and wrote the manuscript.

Competing interests: The authors declare that they have no competing interests.

Data and materials availability: All data necessary to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Material. Additional data related to this paper may be requested from the authors.