Ion transport measurements are widely used as an indirect probe for various properties of confined electrolytes. It is generally assumed that the ion concentration in a nanoscale channel is equal to the ion concentration in the macroscopic reservoirs it connects to, with deviations arising only in the presence of surface charges on the channel walls. Here, we show that this assumption may break down even in a neutral channel due to electrostatic correlations between the ions arising in the regime of interaction confinement, where Coulomb interactions are reinforced due to the presence of the channel walls. We focus on a one-dimensional channel geometry, where an exact evaluation of the electrolyte’s partition function is possible with a transfer operator approach. Our exact solution reveals that in nanometer-scale channels, the ion concentration is generally lower than in reservoirs and depends continuously on the bulk salt concentration, in contrast to the conventional mean-field theory that predicts an abrupt filling transition. We develop a modified mean-field theory taking into account the presence of ion pairs that agrees quantitatively with the exact solution and provides predictions for experimentally relevant observables, such as the ionic conductivity. Our results will guide the interpretation of nanoscale ion transport measurements.

A channel connects two reservoirs filled with a salt solution at concentration cout. What is the salt concentration cin inside the channel? The straightforward answer cin = cout is challenged as soon as the channel’s dimensions are at the nanometer scale.1 A deviation typically occurs because of the presence of a surface charge density Σ on the channel walls, which results in an imbalance of the concentrations cin± of the positive and negative ions. In a cylindrical channel of radius R that is smaller than the electrolyte’s Debye length, under the assumption of local electroneutrality, the concentrations are given by the famous Donnan equilibrium result,2 

cin±=cout2+(2Σ/R)2±2Σ/R.
(1)

Equation (1) is widely used to infer a channel’s surface charge from measurements of its conductivity at different salt concentrations. For sufficiently small surface charges (2Σ/Rcout), Eq. (1) predicts cin = cout even at extreme nanoscales. Importantly, this prediction underlies the method for extracting confined ion mobilities from transport measurements, which has been applied down to 7-Å-wide two-dimensional channels.3 Yet, physically, cin = cout stems from the assumption that the electrolyte solutions, both in the reservoirs and in the channel, behave as ideal gases of non-interacting ions. While such a description is valid in the bulk reservoirs at reasonable salt concentrations,4 it must be challenged in the nanometer-scale channel that is subject to interaction confinement5—a reinforcement of the effective Coulomb interactions between the ions due to the dielectric contrast between the solvent (water) and the channel wall.2,5–13

Due to interaction confinement, ions face a self-energy barrierEs when entering the channel.6,7 It was first noted by Parsegian6 that this should result in ion exclusion: the salt concentration within the channel is then given by an Arrhenius scaling cin=couteEs/kBT under the assumption of non-interacting ions. However, the result becomes more subtle as the confinement-reinforced ionic interactions are taken into account. Within a mean-field description of a spherical nanopore, Dresner14 predicted an abrupt filling transition, where cin was a discontinuous function of cout. Later, Palmeri and co-workers15,16 recovered a similar transition using a three-dimensional model of a cylindrical channel, treated within the variational field theory formalism of Netz and Orland.17 While this approach could be applied to a realistic geometry, it took electrostatic correlations into account only approximately.

An exact treatment of electrostatic correlations is possible upon simplification of the geometry to a purely one-dimensional model, with the channel wall being taken into account by introducing an effective confined Coulomb interaction. The 1D electrolyte can then be mapped onto an Ising or 1D Coulomb-gas-type model; the transfer matrix solution of such models was used, for example, to discuss the capacitance of nanoporous systems.18–20 The lattice models may be taken to the continuum limit, and the resulting path integral solutions have been used to discuss various ion-exchange phase transitions that arise in the presence of fixed discrete charges inside the channel8,21,22 and the ionic Coulomb blockade phenomenon.12 Such models are particularly rich theoretically, as they support a mapping to non-Hermitian quantum mechanics.23 Nevertheless, to our knowledge, the fundamental problem of ion filling in an uncharged channel has not been tackled within this framework.

In this paper, we treat the ion filling problem in the interaction confinement regime using an exactly solvable one-dimensional model. We find that the value of cin is strongly affected by the formation of Bjerrum pairs—pairs of oppositely charged ions—within the channel, which precludes the occurrence of an abrupt filling transition. This is in contrast to the prediction of Palmeri and co-workers15,16 and to the result of the conventional mean-field theory. We then build on our exact results to propose a modified mean-field model that accounts for the relevant physical ingredients and, particularly, for the presence of ion pairs.

The paper is organized as follows: in Sec. II, we present the one-dimensional model and its solution within a path-integral formalism. The reader interested only in the physical outcomes may skip directly to Sec. III, where we discuss the model’s prediction for the ion concentration within the channel, compare it to the mean-field solution, and interpret it in terms of tightly bound Bjerrum pairs. In Sec. IV, we establish a modified mean-field theory, based on the notion of phantom pairs, which reproduces our exact solution. The mean-field theory allows us to determine the number of unpaired ions and produces experimentally relevant predictions for a nanochannel’s ionic conductance. Section V establishes our conclusions.

We consider a cylindrical channel of radius R and length L, connected to macroscopic reservoirs [Fig. 1(a)]. We first assume, for simplicity, that the channel is filled with water that has an isotropic dielectric permittivity ϵw = 80 and that it is embedded in an insulating medium with a much lower permittivity ϵm (for a lipid membrane,6ϵm ∼ 2). The effective Coulomb interaction V(x) between two monovalent ions separated by a distance x on the channel axis can be computed exactly by solving Poisson’s equation.7,11,12 A simple approximate expression can be obtained for xR (Ref. 2),

V(x)e2α2πϵ0ϵwRe|x|/(αR),
(2)

where α is a numerical coefficient that depends on the ratio ϵw/ϵm (α = 6.3 for ϵw/ϵm = 40). The reinforcement of electrostatic interactions compared to the usual e2/4πϵ0ϵwr Coulomb interaction that ions experience in bulk water can be interpreted in terms of image charges within the channel walls [Fig. 1(b)]. Two confined ions interact not only with each other but also with their respective image charges.

FIG. 1.

Ion filling in the interaction confinement regime. (a) Schematic of the ion filling problem: a cylindrical nanochannel (radius R ∼ 1 nm) is connected to macroscopic reservoirs of aqueous electrolyte. The salt concentration inside the channel, cin, may differ from that in the reservoirs, cout. (b) Physics of interaction confinement. When a charged species enters a nanochannel, the dielectric contrast between water (ϵw ∼ 80) and walls (ϵm ∼ 2) constrains the electric field lines to remain within the channel. This process can be interpreted in terms of image charges inside the channel walls and results in an electrostatic self-energy barrier for ions to enter the channel and reinforced interactions between ions.

FIG. 1.

Ion filling in the interaction confinement regime. (a) Schematic of the ion filling problem: a cylindrical nanochannel (radius R ∼ 1 nm) is connected to macroscopic reservoirs of aqueous electrolyte. The salt concentration inside the channel, cin, may differ from that in the reservoirs, cout. (b) Physics of interaction confinement. When a charged species enters a nanochannel, the dielectric contrast between water (ϵw ∼ 80) and walls (ϵm ∼ 2) constrains the electric field lines to remain within the channel. This process can be interpreted in terms of image charges inside the channel walls and results in an electrostatic self-energy barrier for ions to enter the channel and reinforced interactions between ions.

Close modal

We introduce the parameters ξαR and xT ≡ 2πϵ0ϵwR2kBT/e2: both have the dimension of a length. With these notations,

V(x)=kBTξxTe|x|/ξ.
(3)

The effects of ion valence and the anisotropic dielectric response of confined water can be taken into account by adjusting ξ and xT.12 The dielectric response of water is expected to become anisotropic in channels with radius smaller than 5nm:11,24 the component perpendicular to the channel axis is reduced, while the parallel component is largely unaffected. This results in ionic interactions that are even stronger than predicted by Eq. (2), described, in the framework of Eq. (3), by a larger ξ and smaller xT.12 

Formally, the expression in Eq. (2) is valid for any channel radius. Yet, it is only physically relevant if at xR, the interaction is significant compared to kBT, which restricts, in practice, the applicability of Eq. (2) to R ≲ 2 nm. In such extreme 1D confinement, we may neglect the ions’ degrees of freedom perpendicular to the channel axis and assume that they are constrained to move in one dimension. The partition function of such a 1D electrolyte may be computed exactly, as detailed in Sec. II B.

Here, we detail the analytical solution for the partition function of a 1D Coulomb gas-like system that was first introduced in Ref. 12. We set kBT = 1 until the end of Sec. II. We start from a lattice model in order to rigorously establish a path integral description in the continuum limit.

Our computation is inspired by the original solution of the 1D Coulomb gas model by Lenard and Edwards25 and subsequent studies by Démery, Dean, and co-workers,18,20,26,27 and Shklovskii and co-workers.21,22 We consider a one-dimensional lattice with sites 1, …, M as a model for the nanochannel of radius R and length L. Each lattice site i carries a spin Si, which takes the values {0, 1, −1}, corresponding to no ion, a positive ion, or a negative ion occupying the site. We model the surface charge distribution as an extra fixed charge qi added at each lattice site. The spins interact with the Hamiltonian

H({Si})=ξ2xTi,j=1M(Si+qi)(Sj+qj)e|ij|/ξ12xT(S+q)TC(S+q).
(4)

The system is in contact with a particle reservoir at concentration cout. Here, the parameters ξ and xT are dimensionless and expressed in the number of lattice sites.

The grand partition function is given by

Ξ=S1,,SMzi|Si|e12xT(S+q)TC(S+q),
(5)

with z = coutπR2L/M the fugacity. The matrix C can be analytically inverted,

C1=12ξsinh(1/ξ)e1/ξ1000012cosh(1/ξ)100000012cosh(1/ξ)10001e1/ξ.
(6)

Hence, we can carry out a Hubbard–Stratonovich transformation; that is, rewrite the partition function as a Gaussian integral, introducing the integration variable φ,

Ξ=xTM(2π)Mdet(C)S1,,SMzi|Si|dφexT2φTC1φ+i(S+q)Tφ,
(7)

with det(C)=e1/ξ2sinh(1/ξ)ξ(1e2/ξ)M. After performing the sum over the spins, which is now decoupled, we obtain

Ξ=xTM(2π)Mdet(C)dφ1dφMj=1M(1+2zcosφj)j=1MeiqjφjexpxT4ξsinh(1/ξ)j=1M1(φj+1φj)2+2(cosh(1/ξ)1)×j=2M1φj2+(e1/ξ1)(φ12+φM2).
(8)

We now take a continuum limit of the lattice model. We call a the physical lattice spacing and let ξ̃=aξ, x̃T=axT, and z̃=Mz/L. We then let a → 0 and M → ∞ while keeping the physical length of the system L = aM constant. We then drop the tilde sign to lighten the notation and obtain

Ξ=dφ(0)exTφ(0)2/4ξ[dφ]eS[φ]dφ(L)exTφ(L)2/4ξ,
(9)

with

S[φ]=0LdxxT4dφdx2+xT4ξ2φ(x)2iq(x)φ(x)2zcosφ(x)0LL(φ,φ̇).
(10)

q(x) is the one-dimensional density corresponding to the surface charge, and zπR2cout. At this point, ξ and xT have the dimension of length. The path integral measure is defined as

[dφ]=lima0ML=aMj=1MxT4πadφj.
(11)

We now define the propagator P(φ, x|φ0, 0), or simply P(φ, x), as

P(φ,x)=dφ(x)δ(φ(x)φ)[dφ]e0xL(φ,φ̇)×dφ(0)δ(φ(0)φ0).
(12)

Considering an infinitesimal displacement Δx,

P(φ,x)=xT4πΔxd(Δφ)P(φΔφ,xΔx)expxΔxxdxxT4ΔφΔx2+xT4ξ2φ2iq(x)φ2zcosφ.
(13)

Expanding the propagator as P(φ − Δφ, x − Δx) = P(φ, x) − ΔxP/∂x − ΔφP/∂φ + (1/2) (Δφ2)∂2P/∂φ2, and carrying out the Gaussian integrals, we obtain

P(φ,x)=P(φ,x)ΔxPx+O(Δx2)×1ΔxxT4ξ2φ2iq(x)φ2zcosφ+O(Δx2)+ΔxxT2Px2(1+O(Δx)).
(14)

P(φ, x),thus, solves the partial differential equation

Px=1xT2Pφ2+iqφxT4ξ2φ2+2zcosφP,
(15)

with initial condition P(φ, 0) = δ(φφ0), which is the equivalent of a Schrödinger equation for the path integral representation (9). The partition function can, thus, be computed as

Ξ=dφ(L)exTφ2/4ξP(φ,L|f0),
(16)

where P(φ, L|f0) is the solution of (15) with initial condition P(φ,0)=f0(φ)exTφ2/4ξ.

We introduce the Fourier transform of P with respect to φ,

P̃(k,x)=12πdφeikφP(φ,x).
(17)

Then, P̃(k,x)satisfies

P̃x=k2xTP̃qP̃k+xT4ξ22P̃k2+zP̃(k+1,x)+P̃(k1,x).
(18)

From now on, we restrict ourselves to an uncharged channel (q = 0). We then define the operator T such that

[T(P̃)](k)=k2xTP̃+xT4ξ22P̃k2+zP̃(k+1,x)+P̃(k1,x),
(19)

which plays the role of a functional transfer matrix. Recalling Eq. (16), the partition function then reads

Ξ=f0|eLT|f0,
(20)

with f0(k)=eξk2/xT and ⟨f(k)|g(k)⟩ ≡ ∫dkf*(k)g(k).

Now, in the limit L → ∞, we may consider the largest eigenvalue λ of the operator T and the associated eigenfunction χ,

[T(χ)](k)=λχ(k).
(21)

Then, up to an exponentially small correction,

Ξ=|f0|χ|2χ|χeλL.
(22)

Our aim is to compute the salt concentration cin in the nanoscale channel, given a salt concentration cout in the reservoir. At the level of the lattice model, the probability to find, say, a positive ion at position k can be computed by replacing a factor (1 + 2z cos φk) by zeiφk in Eq. (8). In the continuum limit, we obtain the positive (negative) ion linear density at position x by inserting the operator ze (ze) at position x,

πR2cin±(x)=1Ξdφ(0)dφ(x)dφ(L)exTφ(0)2/4ξP(φ(x),x|φ(0),0)×ze±iφ(x)P(φ(L),L|φ(x),x)exTφ(L)2/4ξ.
(23)

Upon Fourier transformation, the insertion of e amounts to a shift by unity. Introducing the operator,

SQ:f(g:kf(kQ)),
(24)

the concentrations are given by

cin±(x)=zπR2f0|exTS±1e(Lx)T|f0Ξ=coutf0|exTS±1e(Lx)T|f0Ξ,
(25)

since z = coutπR2. In the thermodynamic limit, and using Eq. (22) for the partition function, we obtain

cin±=coutχ(k)|χ(k1)χ(k)|χ(k).
(26)

Equation (26) is the main result of our exact computation. In practice, the function χ(k) is determined numerically by the finite-difference integration of Eq. (18).

We now go back to the ion filling problem [Fig. 1(a)] and present first a one-dimensional mean-field solution. Typically, the mean-field solution of an electrolyte problem is obtained by solving the Poisson–Boltzmann equation.28,29 For the conventional Poisson–Boltzmann equation to apply, we would need to consider the full three-dimensional geometry of our problem, and the effective interaction of Eq. (3) would be introduced implicitly through the boundary conditions at the channel walls.14 In order to obtain a mean-field solution directly in the 1D geometry, we need to introduce a modified Poisson’s equation for the electrostatic potential Φ whose Green’s function coincides with Eq. (3),

d2dx21ξ2ϕ=2πR2c+cxT,
(27)

with ϕeΦ/kBT the dimensionless potential. Imposing that the ions follow a Boltzmann distribution (c± = cineϕ, where cin is understood as the average concentration inside the channel), we obtain the analog of the Poisson–Boltzmann equation in our 1D geometry,

d2dx21ξ2ϕ=4πR2cinxTsinhϕ.
(28)

In order to proceed analytically, we make a Debye–Hückel-type approximation and linearize Eq. (28) with respect to ϕ. Then, the potential around an ion placed in the channel at x = 0 is given by

ϕ(x)=ξeffxTe|x|/ξeff,
(29)

with

ξeff2=ξ21+4πR2cinξ2/xT.
(30)

The chemical potential inside the channel is the sum of an ideal gas entropic part and an excess part due to interactions,

μin=μent+μex,
(31)

with

μent=kB TlogcoutΛ3,
(32)

Λ being the de Broglie thermal wavelength of the ions. μex can be obtained via a Debye charging process,30 

μexkBT=01ϕλ(0)dλ,ϕλ(0)=λξ/xT1+4λπR2cinξ2/xT.
(33)

We determine cin by imposing equality of the chemical potentials between the channel and the reservoir,

μout=kBTlogcout Λ3=μin,
(34)

which yields

cin=couteμex/kBT.
(35)

Evaluating analytically the integral in Eq. (33), we obtain an implicit equation for cin. With the notation ĉinπR2cin,

cin=coutexpξ2xT×xT26ξ2ĉin2ξ2132(1+4ĉinξ2/xT)1/2+12(1+4ĉinξ2/xT)3/2.
(36)

In Figs. 2(b) and 2(c), we plot the ratio cin/cout as a function of cout, as obtained by numerically solving Eq. (36). We fix ξ = 7 nm (which corresponds to a channel with R ≈ 1 nm and strong dielectric contrast) and vary xT to set the ionic interaction strength. The interaction strength may be quantified through the self-energy barrier, Es = kBT × ξ/(2xT). The limiting behavior of cin/cout may be understood directly from Eq. (36). When cin is small, Eq. (36) reduces to the Arrhenius scaling cin=couteEs/kB T: this result typically holds for biological ion channels that may contain either 0 or 1 ion at any given time, and the effect of inter-ionic interactions is negligible. When cin is large, we recover cin = cout. Indeed, the excess term in the chemical potential vanishes at high concentrations, which is then dominated by the entropic term. The fact that μex → 0 as cin → ∞ is non-trivial: it can be seen, physically, as resulting from the Coulomb potential of each ion being perfectly screened by the other ions. At small values of Es, Eq. (36) has a single solution for all values of cout, which interpolates smoothly between the two limiting regimes. However, for Es ≳ 5kBT, it has three solutions in a certain range of cout, pointing to a pseudo-first-order phase transition between a low-concentration and a high-concentration phase, similar to the one predicted by Dresner14 and Buyukdagli et al.15 The transition occurs at ĉinxT/ξ2; as per Eq. (30), this corresponds to the concentration where the effect of the screening cloud on an ion’s Coulomb potential becomes significant.

FIG. 2.

Comparing mean-field approximations with the exact Coulomb gas solution. (a) Schematic description of the mean-field approaches. The chemical potential of confined ions is determined by solving the (linear or nonlinear) Poisson–Boltzmann equation around a given ion, interacting with an oppositely charged Debye cloud. (b) Dependence of the channel salt concentration cin on the reservoir salt concentration cout, in a weakly interacting case (R = 1 nm, ξ = 7 nm, xT = 7 nm, and Es = 0.5 kBT). We plot four different predictions for the ratio cin/cout: the exact field-theoretical solution [Eq. (26), blue circles], its low concentration expansion [Eq. (47), black line], the mean-field predictions from solving the full Poisson–Boltzmann equation [Eq. (40), orange curve], or the mean-field predictions from its Debye–Hückel linearization [Eq. (36), yellow line]. The two mean-field predictions are indistinguishable. In all cases, the naive estimate cin = cout is recovered for high enough concentrations. In the dilute limit, the concentration inside the channel is well approximated by the Arrhenius scaling cin=couteEs/kBT. (c) Dependence of the channel salt concentration cin on the reservoir salt concentration cout, in a strongly interacting case (R = 1 nm, ξ = 7 nm, xT = 0.6 nm, and Es = 6 kBT). The color code is the same as in (b). Here, the mean-field predictions strongly deviate from the exact solution, with the Debye–Hückel model predicting an abrupt filling transition. This discrepancy is due to the formation of Bjerrum pairs at intermediate concentrations, as evidenced by the scaling cincout2 in the exact solution.

FIG. 2.

Comparing mean-field approximations with the exact Coulomb gas solution. (a) Schematic description of the mean-field approaches. The chemical potential of confined ions is determined by solving the (linear or nonlinear) Poisson–Boltzmann equation around a given ion, interacting with an oppositely charged Debye cloud. (b) Dependence of the channel salt concentration cin on the reservoir salt concentration cout, in a weakly interacting case (R = 1 nm, ξ = 7 nm, xT = 7 nm, and Es = 0.5 kBT). We plot four different predictions for the ratio cin/cout: the exact field-theoretical solution [Eq. (26), blue circles], its low concentration expansion [Eq. (47), black line], the mean-field predictions from solving the full Poisson–Boltzmann equation [Eq. (40), orange curve], or the mean-field predictions from its Debye–Hückel linearization [Eq. (36), yellow line]. The two mean-field predictions are indistinguishable. In all cases, the naive estimate cin = cout is recovered for high enough concentrations. In the dilute limit, the concentration inside the channel is well approximated by the Arrhenius scaling cin=couteEs/kBT. (c) Dependence of the channel salt concentration cin on the reservoir salt concentration cout, in a strongly interacting case (R = 1 nm, ξ = 7 nm, xT = 0.6 nm, and Es = 6 kBT). The color code is the same as in (b). Here, the mean-field predictions strongly deviate from the exact solution, with the Debye–Hückel model predicting an abrupt filling transition. This discrepancy is due to the formation of Bjerrum pairs at intermediate concentrations, as evidenced by the scaling cincout2 in the exact solution.

Close modal

The physical content of the mean-field solution presented above is similar to the one of Dresner’s solution, based on a linearized Poisson–Boltzmann equation.14 The difference in geometry and the fact that he foregoes the use of the Debye charging process do not seem to play a significant qualitative role. The solution of Buyukdagli et al.15 takes ionic correlations into account to some extent, yet it still involves a Debye–Hückel-type linear equation for the mean-field interaction potential between the ions.

One may ask whether the same phenomenology persists if one does not linearize the Poisson–Boltzmann equation. The full Poisson–Boltzmann equation cannot be solved analytically but supports the following integral form:

dϕdx21ξ2ϕ2=8πR2cinxTcoshϕ1,
(37)

where we have used the fact that ϕ should vanish at x → ∞. For x → 0, the solution of Eq. (37) should reduce to the unscreened potential in Eq. (3) up to an additive constant so that

1xT21ξ2ϕ2(0)=8πR2cinxTcoshϕ(0)1.
(38)

Once again, one may express the excess chemical potential of the confined ions through a Debye charging process,

μexkBT=01ϕλ(0)dλ,λ2xT21ξ2ϕλ2(0)=8πR2λcinxTcoshϕλ(0)1.
(39)

This result is the analog of Eq. (33), with ϕλ(0) now being the solution of an implicit non-linear equation, so that μex must be determined numerically. As before, the concentration inside the channel is then given by

cin=couteμex/kBT.
(40)

The prediction of the full Poisson–Boltzmann equation is shown in Figs. 2(b) and 2(c): we find cin to be a smooth function of cout for all values of parameters, in contrast to the linearized solution. We may not, however, unambiguously conclude that the filling transition is an artifact of linearization since the non-linear solution still involves a mean-field approximation and is not guaranteed to yield the correct result.

Interestingly, the “physically motivated” mean-field solution in Eq. (28) differs from the mean-field limit of our exact solution. It is obtained by taking the saddle-point approximation in the path-integral expression of the partition function [Eq. (9)]. The Euler–Lagrange equation for the minimizer φ(x) of the action S[φ] in Eq. (10) is, upon identifying ϕ = −,

d2dx21ξ2ϕ=4πR2coutxTsinhϕ.
(41)

This is Eq. (28) with cin replaced with cout and corresponds to a first order treatment of interactions. Indeed, if the ions are non-interacting, cin = cout. By solving the mean-field equation, we determine how the ions’ chemical potential is affected by Debye screening, which then results in a value of cin that is different from cout. Within a straightforward interaction expansion procedure, one should determine the effect of screening assuming the zeroth order value for the ion concentration inside the channel, which is cout: this corresponds to Eq. (41). Equation (28) contains an additional self-consistency condition as it assumes the actual value of cin for the ion concentration, which is not known until Eq. (28) is solved. One may draw a loose condensed matter physics analogy, where Eq. (41) resembles the Born approximation for impurity scattering, while Eq. (28) is analogous to the self-consistent Born approximation.31 

We now turn to the exact solution obtained in Sec. II to unambiguously solve the ion filling problem. We determine cin according to Eq. (26) as

cin±=coutχ(k)|χ(k1)χ(k)|χ(k),
(42)

where χ(k) is the highest eigenfunction of the transfer operator in Eq. (19), determined, in practice, by numerical integration. The exact results for cin, with the same parameter values as for the mean-field solution, are shown in Figs. 2(b) and 2(c). When interactions are weak [small values of Es, Fig. 2(b)], the exact and mean-field solutions are in good agreement. Notably, all solutions smoothly interpolate between the bulk scaling cin = cout at a high concentration and the Arrhenius scaling cin=couteEs/kBT at a low concentration. Conversely, in the strongly interacting case [large Es, Fig. 2(c)], the exact result yields a much larger ion concentration than the mean-field solutions for intermediate values of cout. In this intermediate regime, cin remains a smooth function of cout and obeys the scaling cincout2.

Such a scaling is the signature of the formation of tightly bound Bjerrum pairs of positive and negative ions—strongly correlated configurations that are not taken into account by mean-field solutions. Indeed, let us assume that the channel contains an ideal gas of ion pairs at concentration cin. We further assume that in a pair, the distance between the two ions is uniformly distributed in the interval [−xT/2, xT/2] and the binding energy of a pair is kB/xT = 2Es. Then, the grand partition function reads

Ξ=N(zeβEs)2N1N!i=1NLxT/2xT/2dxe2βEs
(43)
=N(z2LxT)NN!=ez2LxT,
(44)

where we recall that z = πR2cout and β ≡ 1/(kBT). Using that

πR2cin=1LlogΞ(βμ)=zLlogΞz,
(45)

we obtain

cin=2z2xTπR2=2πR2xTcout2.
(46)

We recover, indeed, the quadratic scaling.

We may check that the prefactor in Eq. (46) is the correct one by evaluating analytically the expression in Eq. (26) in the low concentration limit zTzxT ≪ 1. An analytical expansion of the function χ(k) in powers of zT was derived in Ref. 12. Substituting it into Eq. (26), we obtain

πR2cin=zeβEs+2zT132zT2eβEs7zT3+O(zT4)+O(e2βEs).
(47)

The first term in the expansion corresponds to cin=couteβEs. At the lowest salt concentrations, forming Bjerrum pairs is too entropically unfavorable and the concentration inside the channel is controlled by the self-energy barrier. However, as the salt concentration increases, there is no abrupt transition to a highly screened concentrated phase inside the channel; instead, the channel is progressively filled by Bjerrum pairs. This corresponds to the quadratic term in the expansion, with the prefactor agreeing, indeed, with Eq. (46). (This justifies a posteriori our choice of [−xT/2, xT/2] as the interval in which a paired-up ion is allowed to move.) The expansion in Eq. (47) reproduces quite well the low-concentration behavior of the exact solution as shown in Figs. 2(b) and 2(c). However, it fails at high concentrations, where it does not recover cin = cout.

Our exact analysis of the ion statistics in a nanoscale channel has revealed that Bjerrum pairs are a crucial ingredient of the filling process. We now develop a modified mean-field theory that accounts for the presence of Bjerrum pairs and compares it to the exact solution.

The traditional mean-field treatment of electrolytes is incapable of taking Bjerrum pairs into account, as it naturally neglects any strong ion–ion correlations—pairing being a fundamentally discrete phenomenon. An idea proposed by Bjerrum to amend the Debye–Hückel theory was to introduce ion pairs as a separate species encapsulating all “strong” ion–ion correlations.32,33 More precisely, any two oppositely charged ions that are closer than some minimum distance can be considered a single neutral entity—a Bjerrum pair. The remaining “free” ions should then only experience weak interactions with each other and can be treated at the mean-field level. Importantly, this last remark justifies the Debye–Hückel linearization, as all non-linear effects are assumed to be hidden in the definition of ion pairs.

As before, we consider that pairs behave like particles of an ideal gas and that their maximum extension is given by xT. Defining cinp the concentration pairs inside the channel, the chemical potential of pairs is given by

μinp=kBTlogcinpΛ62πxTR2,
(48)

where the geometrical factor inside the logarithm accounts for the internal degrees of freedom of a pair. The chemical potential only has an entropic term because the binding energy of the pair exactly compensates for the self-energy of the two separate ions. The chemical equilibrium between free ions and pairs inside the channel can be written as

μin++μin=2μin=μinp,
(49)

where μin+ and μin are the chemical potentials of cations and anions, respectively. We then obtain, using the Debye–Hückel solution for μin [Eqs. (31)(33)],

cinp=2πR2xTcout2,
(50)

which is the result obtained in Sec. III C. The average concentration of free ions cinf is not modified compared to the Debye–Hückel solution and is, therefore, the solution of the self-consistent Eq. (36). One can then compute the total concentration inside the channel as cin=cinf+cinp or, explicitly,

cin=couteμex(cinf)/kBT+2πR2xTcout2.
(51)

In other words, the only impact of pairs in Bjerrum’s computation is to add a quadratic term 2πR2xTcout2 to the Debye–Hückel result, matching with the expansion (47) of the exact solution up to order 2 in the bulk concentration. We compare the two predictions in Fig. 3(b). The Debye–Hückel–Bjerrum solution is found to match the exact one quite well at low and intermediate concentrations. This result is, however, unphysical for cout ≳ 1/πR2xT: cin is found to grow much faster than the bulk concentration. One solution would be to consider higher-order terms in the mean-field treatment through the inclusion of triplets, quadruplets, etc., of ions and all possible interactions between these entities. Truncating the sum at any finite order, however, would not yield a valid solution in the entire range of concentrations, nor is it guaranteed to converge to the exact solution. This approach is also unsatisfactory as it would not yield a closed-form expression for cin and would not allow for a qualitative understanding of the underlying physics.

FIG. 3.

Pair-enhanced mean-field theory. (a) Treatment of ion pairing in mean-field approaches. Top panel: Mean-field theories inevitably underestimate ion–ion correlations. To circumvent this problem, two ions that are distant by less than xT are considered to form an ion pair, which is treated as a separate chemical species. Bottom panel: schematic representation of ion distribution around a fixed positive ion. The distribution is very peaked close to the central ion, due to the formation of an ion pair, and then relaxes smoothly to the mean value cin. (b) Evolution of channel concentration cin as a function of reservoir concentration cout, in a strongly interacting case (R = 1 nm, ξ = 7 nm, xT = 0.6 nm, and Es = 6 kBT). We plot the ratio cin/cout obtained from three different models taking Bjerrum pairs into account: the exact field-theoretical solution [Eq. (26), blue circles], the Debye–Hückel–Bjerrum mean-field theory [Eq. (51), red line], and our modified mean-field theory based on the notion of phantom pairs [Eq. (55), orange line], which reproduces the exact solution quantitatively for all values of parameters. At high concentrations, the Debye–Hückel–Bjerrum prediction fails due to the uncontrolled proliferation of Bjerrum pairs. (c) Pairs and phantom pairs. Top panel: In the channel, oppositely charged ions may form tightly bound Bjerrum pairs due to strong Coulomb interactions. Bottom panel: We equilibrate the channel with a fictitious 1D reservoir of non-interacting ions. Particularly at high concentrations, two oppositely charged ions may find themselves within a distance xT out of statistical chance: this is a phantom pair.

FIG. 3.

Pair-enhanced mean-field theory. (a) Treatment of ion pairing in mean-field approaches. Top panel: Mean-field theories inevitably underestimate ion–ion correlations. To circumvent this problem, two ions that are distant by less than xT are considered to form an ion pair, which is treated as a separate chemical species. Bottom panel: schematic representation of ion distribution around a fixed positive ion. The distribution is very peaked close to the central ion, due to the formation of an ion pair, and then relaxes smoothly to the mean value cin. (b) Evolution of channel concentration cin as a function of reservoir concentration cout, in a strongly interacting case (R = 1 nm, ξ = 7 nm, xT = 0.6 nm, and Es = 6 kBT). We plot the ratio cin/cout obtained from three different models taking Bjerrum pairs into account: the exact field-theoretical solution [Eq. (26), blue circles], the Debye–Hückel–Bjerrum mean-field theory [Eq. (51), red line], and our modified mean-field theory based on the notion of phantom pairs [Eq. (55), orange line], which reproduces the exact solution quantitatively for all values of parameters. At high concentrations, the Debye–Hückel–Bjerrum prediction fails due to the uncontrolled proliferation of Bjerrum pairs. (c) Pairs and phantom pairs. Top panel: In the channel, oppositely charged ions may form tightly bound Bjerrum pairs due to strong Coulomb interactions. Bottom panel: We equilibrate the channel with a fictitious 1D reservoir of non-interacting ions. Particularly at high concentrations, two oppositely charged ions may find themselves within a distance xT out of statistical chance: this is a phantom pair.

Close modal

Instead, we develop a different method that, through physics-driven arguments, prevents the divergence of cin at high bulk concentrations and quantitatively reproduces the exact solution.

Equation (51) overestimates the number of Bjerrum pairs in the channel because it fails to account for the presence of Bjerrum pairs in the reservoir. The electrolyte in the reservoir is treated as an ideal gas: the ions are non-interacting and they cannot form actual, tightly bound pairs. Nevertheless, we have defined any two oppositely charged ions that find themselves in a cylinder of radius R and length xT to be a separate chemical species. Such configurations may arise in the reservoir simply out of statistical chance: we dub them phantom pairs. For our mean-field theory to be consistent, these phantom pairs need to be taken into account.

We incorporate phantom pairs in our model using a thought experiment: let us assume that the nanochannel does not directly connect to the reservoir but is, instead, in contact with another 1D channel of similar geometry, where interactions are turned off. This intermediate 1D reservoir is, in turn, connected to the real, 3D reservoir. Let coutp (resp. coutf) be the concentration of the aforementioned phantom pairs (resp. free ions) in the 1D reservoir. The chemical equilibrium between phantom pairs and free ions imposes

coutp=2πR2xT(coutf)2.
(52)

In addition, one has coutf+coutp=cout: since the interactions are negligible, the total concentration of ions in the intermediate reservoir should be equal to the salt concentration in the real reservoir. Imposing this condition yields

coutf=1+8πcoutxTR214xTπR2.
(53)

We use this result to control the proliferation of pairs in the channel: we now equilibrate the free ions inside the nanochannel with only the free ions in the 1D reservoir,

cinf=coutfeμex(cinf)/kBT,
(54)

which corresponds to Eq. (35) with cout replaced by coutf. Equation (54) is again a self-consistent equation; this time on the concentration of free ions cinf that must be solved numerically. Finally, equilibrating pairs with free ions inside the channel (or, equivalently, pairs inside with pairs outside), we obtain

cin=cinf+2πR2xT(coutf)2,
(55)

where the second term corresponds again to Bjerrum pairs. Equations (53)(55) constitute the main result of our modified mean-field theory. Note that μex may be determined at the Debye–Hückel level [Eq. (33)] or by solving the full Poisson–Boltzmann equation [Eq. (39)]. In what follows, we will only discuss the latter, as it offers greater accuracy; however, the Debye–Hückel prediction provides reasonable results even in the case of strong interactions and yields a convenient analytical expression for μex as a function of cinf.

The prediction of our phantom pair Poisson–Boltzmann model is compared to the exact solution (26) in Fig. 3(b). The two solutions are found to be in near perfect agreement for all values of the parameters, even in the strong coupling limit EskBT. This result shows the importance of taking into account phantom pairs in the reservoir to solve the ion filling problem.

In Secs. II C and II D, we use our modified mean-field model to predict the conductance of a nanochannel, first in the case of a neutral channel and then in the presence of a surface charge.

A strength of our modified mean-field model is that it offers insights into the physical properties of the confined system beyond the value of the ionic concentration. In particular, the decomposition of the electrolyte into free ions and bound pairs allows us to estimate the channel’s conductance. Tightly bound Bjerrum pairs are electrically neutral so that they do not contribute to the ionic current to the first order in an applied electric field: it would then be straightforward to assume that the channel’s conductance is proportional to the concentration of free ions. However, the reasoning needs to be more subtle, since the channel, in the same way as the reservoir, may contain non-interacting phantom pairs. Indeed, we have decomposed the confined electrolyte into tightly bound pairs that have no ionic atmosphere and free ions that are dressed by a Debye screening cloud. As the concentration increases, the interaction between dressed ions becomes weak and two of them may find themselves within a distance xT without actually binding. Such a phantom pair is expected to still contribute to the conductance. The concentration of phantom pairs in the channel is obtained by imposing their chemical equilibrium with the free ions treated as an ideal gas. Thus, we estimate the channel’s conductance as

G=2e2DkBTπR2Lcinf+2xTπR2(cinf)2,
(56)

where D is the diffusion coefficient of ions and the second term corresponds to the contribution of phantom pairs. In Fig. 4(a), we compare this result to the Ohm’s law prediction where pairs are neglected and one assumes cin = cout. Ohm’s law is found to greatly overestimate the conductance at low concentrations. In the dilute limit, we, instead, recover the Arrhenius scaling, where one assumes cin=couteEs/kBT.

FIG. 4.

Channel conductance in the pair-enhanced mean-field model. (a) Conductance of a nanochannel (R = 1 nm, ξ = 7 nm, xT = 0.7 nm, and Es = 10 kBT) as a function of the reservoir concentration. The red line corresponds to the prediction of the phantom pair mean-field model [Eq. (56)] for T = 300 K, D = 10−9 m2/s, and L = 100 nm. The Ohm’s law bulk prediction (cin = cout, blue line) and the Arrhenius model (cin=couteEs/kBT, yellow line) are also represented for comparison. (b) Conductance of a nanochannel with a weak surface charge Σ = 10−3 C/m2. We represented the predictions of the conventional Donnan equilibrium [Eq. (1), blue line] and of the phantom pair mean-field theory [Eqs. (56) and (59), red line]. Because interaction confinement results in a lower ion concentration in the channel, the usual formula Σ ∼ Rc*/2, where c* is the reservoir concentration for which the conductance levels off, overestimates the surface charge by one order of magnitude, as indicated on the plot.

FIG. 4.

Channel conductance in the pair-enhanced mean-field model. (a) Conductance of a nanochannel (R = 1 nm, ξ = 7 nm, xT = 0.7 nm, and Es = 10 kBT) as a function of the reservoir concentration. The red line corresponds to the prediction of the phantom pair mean-field model [Eq. (56)] for T = 300 K, D = 10−9 m2/s, and L = 100 nm. The Ohm’s law bulk prediction (cin = cout, blue line) and the Arrhenius model (cin=couteEs/kBT, yellow line) are also represented for comparison. (b) Conductance of a nanochannel with a weak surface charge Σ = 10−3 C/m2. We represented the predictions of the conventional Donnan equilibrium [Eq. (1), blue line] and of the phantom pair mean-field theory [Eqs. (56) and (59), red line]. Because interaction confinement results in a lower ion concentration in the channel, the usual formula Σ ∼ Rc*/2, where c* is the reservoir concentration for which the conductance levels off, overestimates the surface charge by one order of magnitude, as indicated on the plot.

Close modal

Finally, we stress that Eq. (56) only accounts for the electrophoresis of free ions and is, therefore, only valid in the limit of weak external electric fields. Stronger voltage drops will result in the breaking of ion pairs, causing a conductivity increase in a process known as the second Wien effect. This phenomenon is described in Refs. 12 and 13 and has been used to create solid-state voltage-gated nanochannels.34 Furthermore, stronger electric fields may result in an unconventional transport mechanism involving the breakdown of the water-ion chain inside the channel, in apparent violation of the Nernst–Einstein relation.35 

Until now, we have restricted ourselves to channels with uncharged walls. However, in most experimentally relevant situations, the channel walls bear a surface charge density Σ, which strongly impacts nanofluidic transport. While introducing a surface charge is tedious within the exact framework, we may readily assess the effect of a surface charge in the interaction confinement regime using our pair-enhanced mean-field theory.

In the limit where the channel’s radius is smaller than the Debye length, we assume that the presence of the surface charge amounts to a homogeneous Donnan potential drop VD inside the channel, which we do not need to determine explicitly. Then, the chemical potential of ions inside the channel reads

μin±=μex±eVD+kBTlogcin±Λ3.
(57)

Note that the concentration in free anions cin and cations cin+ are now distinct so that μex is defined as a function of the average free ion concentration cinf=(cin++cin)/2. In a channel of finite length, the ion concentrations do not necessarily obey local electroneutrality.36 It may, indeed, be more favorable for counterions to form a cloud at the channel mouths instead of entering the channel. At low salt concentrations, this can occur even in very long isolated channels, but the effect is mitigated in a membrane configuration with many parallel channels.37 Such non-neutrality effects are, however, beyond the scope of this work, and we assume in the following that the channel is sufficiently long for local electroneutrality to hold,

cin+cin+2Σ/R=0.
(58)

Imposing chemical equilibrium with the reservoir, we obtain a modified version of the Donnan result [Eq. (1)]:

cin=cinf+cinp,cinf=coutfeβμex(cinf)2+2ΣR2,cinp=2πR2xT(coutf)2,
(59)

with coutf given by Eq. (53).

One can again obtain the channel’s conductance through Eq. (56), which we compare to the Donnan/Ohm’s law result in Fig. 4(b). Importantly, the Donnan result predicts that conductance becomes independent of concentration for cout ∼ 2Σ/R [see Eq. (1)]. In practice, this result is commonly used to estimate experimentally the surface charge as Σ ∼ Rc*/2, where c* is the reservoir concentration for which the conductance levels off. In contrast, in the interaction confinement regime, we predict that the transition occurs, instead, at cinf2Σ/R—corresponding to a higher reservoir concentration, due to the self-energy barrier. In this case, Donnan’s prediction overestimates the surface charge by typically one order of magnitude, as shown in Fig. 4(b).

Finally, let us stress that we considered here a charge homogeneously distributed along the channel’s surface. This assumption is relevant in the case of conducting wall materials, such as systems where the charge is imposed via a gating electrode connected to the channel walls. This situation may, however, be different in experimentally available devices, where the surface charge generally consists in localized charged groups and defects on the channel walls. In this case, the physics become more involved as ions may form bound pairs with the fixed surface charges. Some of these physics have been revealed by the exact computations of Shklovskii and co-workers;8,21 a technically simpler approach to these physics using our pair-enhanced mean-field theory would be possible but extends beyond the scope of the present work.

We have determined the salt concentration inside a nanometric channel connected to reservoirs filled with electrolytes. In the case of a fully 1D geometry, corresponding to a nanotube of radius R ∼ 1 nm, we developed an exact field-theoretical solution that allowed us to compute channel concentration cin as a function of the reservoir concentration cout. This solution clears up the ambiguities of pre-existing mean-field theories and contradicts the naive expectation that cin = cout. In particular, the concentration inside the nanochannel is found to be always lower than in the bulk, as the confinement of electrostatic interactions creates an energy barrier for ions to enter the channel.

Yet, we found that cin is, in fact, higher than the prediction of the mean-field Debye–Hückel theory, as ion pairing counterbalances, to some extent, the energy cost of interaction confinement. Such strong ion–ion correlations cannot be directly accounted for in a mean-field theory, and the filling transition that emerges in the Debye–Hückel theory appears to be an artifact of linearization. To overcome this issue, one can add Bjerrum pairs as a separate chemical species within the Debye–Hückel model. Carefully accounting for the statistical formation of unbound phantom pairs, we obtain a modified mean-field theory that reproduces the result of the exact computation with a near-perfect accuracy and that can be extended to account for a non-zero surface charge on the channel wall.

Despite the concurring results, the two original formalisms developed in this work serve distinct purposes. The field-theoretical solution plays the role of a touchstone model owing to its exact treatment of all many-body interactions. Modeling electrolytes is a notoriously hard problem in statistical physics, and simplified models often lack a reference solution for benchmarking their approximations. This difficulty is lifted in the 1D geometry: thanks to the existence of the exact solution, we have been able to build a quantitatively precise mean-field model, adding step-by-step the qualitative ingredients necessary to reproduce the exact result.

Moreover, the field theory formalism gives access to the entire statistics of the system, including finite-size effects that elude any mean-field treatment. The latter are expected to be relevant in many experimental situations, as a substantial amount of current work focuses on short pores, where the length of the channel is comparable to its radius. For instance, one can expect shorter channels to deviate from electroneutrality36—something entirely impossible in the limit of infinitely long channels.

On the other hand, our modified mean-field formalism has the advantage of mathematical simplicity, allowing for convenient physical interpretations. The simple distinction between free ions and Bjerrum pairs can be used to straightforwardly estimate the channel’s conductance. The influence of ion–ion correlations on the conductivity is of particular importance as conductance measurements underpin many nanofluidic experiments. In contrast, the exact solution does not provide any such insight into transport properties, as it is limited to thermal equilibrium.

Furthermore, the mean-field model may easily be adapted to other geometries, whereas an exact treatment is only possible in the strictly 1D case. Extensions of our results to 2D nanochannels would be of significant interest. In particular, 2D nanochannels can be made out of various materials with different electronic properties, which directly impact the confined ionic interactions.5 Therefore, 2D nanochannels could serve as a platform for exploring the impact of wall metallicity on the ion filling problem.

Both our exact and mean-field solutions can be expected to fail at very high concentrations. Indeed, our work relies on a simplified picture of electrolytes, where all steric effects are discarded. We considered point-like ions with no short-distance repulsion; therefore, no effects such as saturation or layering can be accounted for. Similarly, we neglected any interaction with the solvent—for example, we did not consider the decrement in relative permittivity at high salt concentrations.38 However, since all electrostatic interactions are screened in the limit of high concentrations, such considerations should not impact the conclusions of the present work: particularly, we would still expect that cin = cout at high concentrations.

Finally, let us briefly recall our results for the ion filling problem. In channels larger than a few nanometers, the conventional mean-field picture is valid so that in the absence of any surface charge, the salt concentration inside the channel equals that of the reservoirs: cin = cout. For nanometer-scale confinement and low concentrations, interaction confinement amounts to a finite energy barrier for ions to enter the channel: cin=couteEs/kBT. As concentration increases, more ions are able to overcome the barrier by forming Bjerrum pairs, neutralizing the electrostatic cost of confinement, at the price of entropy: cincout2. Only at high concentrations can one recover the intuitive estimate cin = cout, as intense screening cancels out all electrostatic interactions. Overall, the interaction confinement has a significant impact on the properties of nanofluidic systems and the assumption that cin = cout should be questioned any time the system’s size reaches the nanometer scale.

N.K. acknowledges the support from a Humboldt fellowship. L.B. acknowledges the funding from the EU H2020 Framework Programme/ERC Advanced Grant Agreement No. 785911-Shadoks. The Flatiron Institute is a division of the Simons Foundation.

The authors have no conflicts to disclose.

Paul Robin: Conceptualization (equal); Investigation (equal); Writing – original draft (equal); Writing – review & editing (equal). Adrien Delahais: Investigation (supporting). Lydéric Bocquet: Conceptualization (equal); Writing – review & editing (equal). Nikita Kavokine: Conceptualization (equal); Investigation (equal); Writing – original draft (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.

1.
R. B.
Schoch
,
J.
Han
, and
P.
Renaud
, “
Transport phenomena in nanofluidics
,”
Rev. Mod. Phys.
80
,
839
883
(
2008
).
2.
N.
Kavokine
,
R. R.
Netz
, and
L.
Bocquet
, “
Fluids at the nanoscale: From continuum to subcontinuum transport
,”
Annu. Rev. Fluid Mech.
53
,
377
410
(
2021
).
3.
A.
Esfandiar
,
B.
Radha
,
F. C.
Wang
,
Q.
Yang
,
S.
Hu
,
S.
Garaj
,
R. R.
Nair
,
A. K.
Geim
, and
K.
Gopinadhan
, “
Size effect in ion transport through angstrom-scale slits
,”
Science
358
,
511
513
(
2017
).
4.
Y.
Avni
,
R. M.
Adar
,
D.
Andelman
, and
H.
Orland
, “
Conductivity of concentrated electrolytes
,”
Phys. Rev. Lett.
128
,
098002
(
2022
).
5.
N.
Kavokine
,
P.
Robin
, and
L.
Bocquet
, “
Interaction confinement and electronic screening in two-dimensional nanofluidic channels
,”
J. Chem. Phys.
157
,
114703
(
2022
).
6.
A.
Parsegian
, “
Energy of an ion crossing a low dielectric membrane: Solutions to four relevant electrostatic problems
,”
Nature
221
,
844
846
(
1969
).
7.
S.
Teber
, “
Translocation energy of ions in nano-channels of cell membranes
,”
J. Stat. Mech.: Theory Exp.
2005
,
P07001
.
8.
J.
Zhang
,
A.
Kamenev
, and
B. I.
Shklovskii
, “
Conductance of ion channels and nanopores with charged walls: A toy model
,”
Phys. Rev. Lett.
95
,
148101
(
2005
).
9.
Y.
Levin
, “
Electrostatics of ions inside the nanopores and trans-membrane channels
,”
Europhys. Lett.
76
,
163
169
(
2006
).
10.
S.
Kondrat
and
A.
Kornyshev
, “
Superionic state in double-layer capacitors with nanoporous electrodes
,”
J. Phys.: Condens. Matter
23
,
022201
(
2011
).
11.
P.
Loche
,
C.
Ayaz
,
A.
Schlaich
,
Y.
Uematsu
, and
R. R.
Netz
, “
Giant axial dielectric response in water-fille nanotubes an effective electrostatic ion–ion interactions from a tensorial dielectric model
,”
J. Phys. Chem. B
123
,
10850
10857
(
2019
).
12.
N.
Kavokine
,
S.
Marbach
,
A.
Siria
, and
L.
Bocquet
, “
Ionic Coulomb blockade as a fractional Wien effect
,”
Nat. Nanotechnol.
14
,
573
578
(
2019
).
13.
P.
Robin
,
N.
Kavokine
, and
L.
Bocquet
, “
Modeling of emergent memory and voltage spiking in ionic transport through angstrom-scale slits
,”
Science
373
,
687
691
(
2021
).
14.
L.
Dresner
, “
Ion exclusion from neutral and slightly charged pores
,”
Desalination
15
,
39
57
(
1974
).
15.
S.
Buyukdagli
,
M.
Manghi
, and
J.
Palmeri
, “
Ionic capillary evaporation in weakly charged nanopores
,”
Phys. Rev. Lett.
105
,
158103
(
2010
).
16.
S.
Buyukdagli
,
M.
Manghi
, and
J.
Palmeri
, “
Variational approach for electrolyte solutions: From dielectric interfaces to charged nanopores
,”
Phys. Rev. E
81
,
041601
(
2010
).
17.
R. R.
Netz
and
H.
Orland
, “
Variational charge renormalization in charged systems
,”
Eur. Phys. J. E
11
,
301
311
(
2003
).
18.
V.
Démery
,
D. S.
Dean
,
T. C.
Hammant
,
R. R.
Horgan
, and
R.
Podgornik
, “
The one-dimensional Coulomb lattice fluid capacitor
,”
J. Chem. Phys.
137
,
064901
(
2012
).
19.
A. A.
Lee
,
S.
Kondrat
, and
A. A.
Kornyshev
, “
Single-file charge storage in conducting nanopores
,”
Phys. Rev. Lett.
113
,
048701
(
2014
).
20.
V.
Démery
,
R.
Monsarrat
,
D. S.
Dean
, and
R.
Podgornik
, “
Phase diagram of a bulk 1d lattice Coulomb gas
,”
Europhys. Lett.
113
,
18008
(
2016
).
21.
J.
Zhang
,
A.
Kamenev
, and
B. I.
Shklovskii
, “
Ion exchange phase transitions in water-filled channels with charged walls
,”
Phys. Rev. E
73
,
051205
(
2006
).
22.
A.
Kamenev
,
J.
Zhang
,
A. I.
Larkin
, and
B. I.
Shklovskii
, “
Transport in one-dimensional Coulomb gases: From ion channels to nanopores
,”
Phys. A
359
,
129
161
(
2006
).
23.
T.
Gulden
and
A.
Kamenev
, “
Dynamics of ion channels via non-Hermitian quantum mechanics
,”
Entropy
23
,
125
(
2021
).
24.
L.
Fumagalli
,
A.
Esfandiar
,
R.
Fabregas
,
S.
Hu
,
P.
Ares
,
A.
Janardanan
,
Q.
Yang
,
B.
Radha
,
T.
Taniguchi
,
K.
Watanabe
,
G.
Gomila
,
K. S.
Novoselov
, and
A. K.
Geim
, “
Anomalously low dielectric constant of confined water
,”
Science
360
,
1339
1342
(
2018
).
25.
S. F.
Edwards
and
A.
Lenard
, “
Exact statistical mechanics of a one-dimensional system with Coulomb forces. II. The method of functional integration
,”
J. Math. Phys.
3
,
778
792
(
1962
).
26.
V.
Démery
,
D. S.
Dean
,
T. C.
Hammant
,
R. R.
Horgan
, and
R.
Podgornik
, “
Overscreening in a 1D lattice Coulomb gas model of ionic liquids
,”
Europhys. Lett.
97
,
28004
(
2012
).
27.
D. S.
Dean
,
R. R.
Horgan
, and
D.
Sentenac
, “
Boundary effects in the one-dimensional Coulomb gas
,”
J. Stat. Phys.
90
,
899
926
(
1998
).
28.
D.
Andelman
, “
Electrostatic properties of membranes: The Poisson–Boltzmann theory
,”
Handbook of Biological Physics
(
Elsevier
,
1995
), Vol. 1, pp.
603
642
.
29.
C.
Herrero
and
L.
Joly
, “
Poisson–Boltzmann formulary
,” arXiv:2105.00720 (
2021
).
30.
D.
Frenkel
and
B.
Smit
,
Understanding Molecular Simulation
(
Academic Press
,
2002
), Chap. 7.
31.
H.
Bruus
and
K.
Flensberg
,
Many-Body Quantum Theory in Condensed Matter Physics
(
Oxford University Press
,
2016
), Chap. 12.
32.
Y.
Levin
, “
Electrostatic correlations: From plasma to biology
,”
Rep. Prog. Phys.
65
,
1577
(
2002
).
33.
V.
Neklyudov
and
V.
Freger
, “
Putting together the puzzle of ion transfer in single-digit carbon nanotubes: Mean field meets ab initio
,”
Nanoscale
14
,
8677
8690
(
2022
).
34.
P.
Robin
,
T.
Emmerich
,
A.
Ismail
,
A.
Niguès
,
Y.
You
,
G.-H.
Nam
,
A.
Keerthi
,
A.
Siria
,
A.
Geim
,
B.
Radha
, and
L.
Bocquet
, “
Long-term memory and synapse-like dynamics of ionic carriers in two-dimensional nanofluidic channels
,”
Science
379
,
161
(
2023
).
35.
Z.
Li
,
R. P.
Misra
,
Y.
Li
,
Y.-C.
Yao
,
S.
Zhao
,
Y.
Zhang
,
Y.
Chen
,
D.
Blankschtein
, and
A.
Noy
, “
Breakdown of the Nernst–Einstein relation in carbon nanotube porins
,”
Nat. Nanotechnol.
18
,
177
(
2022
).
36.
A.
Levy
,
J. P.
de Souza
, and
M. Z.
Bazant
, “
Breakdown of electroneutrality in nanopores
,”
J. Colloid Interface Sci.
579
,
162
176
(
2020
).
37.
J. P.
de Souza
,
A.
Levy
, and
M. Z.
Bazant
, “
Electroneutrality breakdown in nanopore arrays
,”
Phys. Rev. E
104
,
044803
(
2021
).
38.
A.
Levy
,
D.
Andelman
, and
H.
Orland
, “
Dipolar Poisson–Boltzmann approach to ionic solutions: A mean field and loop expansion analysis
,”
J. Chem. Phys.
139
,
164909
(
2013
).