Electroosmosis in a finite cylindrical pore: simple models of end effects

TL;DR

This paper develops simple models to analyze electroosmotic flow in finite cylindrical pores, accounting for end effects and charge cloud spilling, with validation against numerical computations across different regimes.

Contribution

It introduces an approximate analysis that combines limiting cases to effectively model electroosmosis in finite pores, including end effects and charge spill-out.

Findings

Full numerical results agree with theory for h>>a and h=0 cases.

Approximate analysis captures main features for small Debye length.

Charge spill-out reduces electroosmotic flow when Debye length is large.

Abstract

A theoretical model of electroosmosis through a circular pore of radius $a$ that traverses a membrane of thickness $h$ is investigated. Both the cylindrical surface of the pore and the outer surfaces of the membrane are charged. When $h\gg a$ end effects are negligible: results of full numerical computations of electroosmosis in an infinite pore agree with theory. When $h=0$, end effects dominate, and computations again agree with analysis. For intermediate values of $h/a$, an approximate analysis that combines these two limiting cases captures the main features of computational results when the Debye length $\kappa^{-1}$ is small compared with the pore radius $a$. However, the approximate analysis fails when $\kappa^{-1}\gg a$, when the charge cloud due to the charged cylindrical walls of the pore spills out of the ends of the pore, and the electroosmotic flow is reduced. When this spilling out is included in the analysis, agreement with computation is restored.

Tables

Table 1. Table 1: The charge lost from the ends of a charged pore when the membrane charge density σ m = 0 subscript 𝜎 𝑚 0 \sigma_{m}=0 , in terms of an equivalent pore length h lost subscript ℎ lost h_{\text{lost}} ( 64 ), and the charge gained inside an uncharged pore ( σ c = 0 subscript 𝜎 𝑐 0 \sigma_{c}=0 ) from the charge cloud adjacent to the charged membrane surfaces, in terms of an equivalent pore length h gained subscript ℎ gained h_{\text{gained}} ( 77 ). a ​ κ = 0.1 𝑎 𝜅 0.1 a\kappa=0.1 .

$h/a$	$h\kappa$	$h_{\text{lost}}/a$		$h_{\text{gained}}/a$
		theory (64)	numerical	theory (77)	numerical
10.0	1.0	8.9186	8.9249	0.4081	0.4119
1.0	0.1	0.9953	0.9954	0.0455	0.0459
0.1	0.01	0.1000	0.1000	0.0046	0.0046

Equations

$$z=a\sinh \xi \cos \eta ,r=a\cosh \xi \sin \eta ,$$

$${\nabla }^2{\varphi }_0={\kappa }^2{\varphi }_0,$$

$${\rho }_0=-ϵ{\kappa }^2{\varphi }_0.$$

$${\nabla }^2\chi =0,$$

$$\mathbf{n}.\nabla \chi =0$$

$$\chi =\frac{{\varphi }_1}{2}\left[1-\frac{2}{\pi }{\tan }^{-1}\left(\frac{1}{\sinh \xi }\right)\right]={\tilde{\chi }}_m(r,z){\varphi }_1.$$

$${\tilde{\chi }}_m=0,z=0,r<a,h=0.$$

$$I_m=-\frac{{\varphi }_1}{R_m},R_m=\frac{1}{2a\Sigma }.$$

$$\chi ={\tilde{\chi }}_c{\varphi }_2=\frac{z}{h}{\varphi }_2,r<a,|z|<h/2,$$

$$\chi =\frac{{\varphi }_2}{2}+({\varphi }_1-{\varphi }_2){\tilde{\chi }}_m(r,z-h/2),$$

$$I_c=-\frac{{\varphi }_2}{R_c},R_c=\frac{h}{\pi a^2\Sigma },$$

$$I_m=-\frac{{\varphi }_1-{\varphi }_2}{R_m}.$$

$${\varphi }_2=\frac{R_c{\varphi }_1}{R_m+R_c}.$$

$${\varphi }_0={\zeta }_c\frac{I_0(\kappa r)}{I_0(\kappa a)}=\frac{{\sigma }_c}{ϵ\kappa }\frac{I_0(\kappa r)}{I_1(\kappa a)}.$$

$$u=\frac{ϵ{\varphi }_2}{\mu h}({\zeta }_c-{\varphi }_0),$$

$$Q_{ce}=\frac{2\pi {\sigma }_c{a}^3}{\mu h}\left[\frac{1}{2\kappa a}\frac{I_0(\kappa a)}{I_1(\kappa a)}-\frac{1}{(\kappa a{)}^2}\right]{\varphi }_2=H_c{\varphi }_2,$$

$$H_{c}&\sim\frac{\pi\sigma_{c}a^{2}}{\mu h\kappa}=\frac{\pi a^{2}\zeta_{c}}{\mu h\epsilon},\hskip 20.0pta\kappa\gg 1,\\ \sim\frac{\pi\sigma_{c}a^{3}}{4\mu h},\hskip 20.0pta\kappa\ll 1.{}$$

$$\mathbf{u}=p_1\mathbf{G},$$

$$Q=-{\int }_V{\rho }_0\mathbf{G}.\nabla \chi \text{d}V,$$

$$\mathbf{u}=p_1{\mathbf{G}}^m.$$

$${\rho }_0={\sigma }_m{\kappa }^{2}a\left[{\int }_0^{\infty }\frac{J_1(as)J_0(rs)}{({\kappa }^2+s^2{)}^{1/2}}{e}^{-({\kappa }^2+s^2{)}^{1/2}z}\text{d}s-\frac{e^{-\kappa z}}{\kappa a}\right],$$

$$Q_{me}=H_m{\varphi }_1,$$

$$H_m\sim a\kappa H_0,a\kappa \ll 1,$$

$$H_0=\frac{a^2{\sigma }_m}{3\mu }.$$

$$Q_{me}\approx H_m({\varphi }_1-{\varphi }_2),$$

$$Q_{ce}\approx H_c{\varphi }_2.$$

$$p=0,z=0,r<a,h=0.$$

$$Q_{mh}=G_m{p}_1,G_m=-\frac{a^3}{3\mu }.$$

$$p=\frac{p_2}{h}z,r<a,|z|<h/2,$$

$$\mathbf{u}=p_2{\mathbf{G}}^c,$$

Full text

\altaffiliation

Department of Engineering Sciences and Applied Mathematics, Northwestern University, 2145 Sheridan Rd, Evanston, IL 60208, USA

{tocentry}

1 1. Introduction

Electroosmosis in a circular cylindrical pore of finite length $h$ differs from that in an infinitely long pore due to end effects. If the cylinder length $h=0$, the pore consists of a hole in a charged membrane of zero thickness, and electroosmosis can be considered to be entirely due to end effects. This case was considered by us previously 1. When the cylindrical pore is infinitely long, end effects are negligible, and the computation of the electroosmotic volumetric flow rate $Q$ is straightforward2, 3. Here we are interested in intermediate values of $h$. Similar results are also available in planar infinitely long channels for arbitrary Debye length and wall charge 4, 5. A related problem of interest is the flow generated by an electric field along a flat surface with a step change in surface charge density 6, 7.

Full numerical computation of the Poisson-Nernst-Planck (PNP) equations for ionic motion is of course possible, and some typical results were reported by Mao et al. 1. Such numerical computations however do not identify the mechanisms underlying the qualitative features of the physical system. Here we discuss how simple models, based upon continuity of electric current and volumetric flow rate, can be combined in order to estimate end effects for pore lengths $h>0$. We assume that the zeta potential on the surface of the membrane is small, so that the Poisson-Boltzmann equation governing the equilibrium charge cloud can be linearized, and the electroosmotic velocity can be determined by an analysis equivalent to that of Henry 14 for electrophoresis, i.e. fluid motion is generated by the effect of the applied electric field acting on the equilibrium charge cloud (which is not deformed either by the applied electric field or by fluid motion). In this limit the electroosmotic volumetric flow rate $Q$ through the hole in the membrane can be determined by means of the reciprocal theorem 1.

1 shows the axisymmetric geometry that we are considering. The cylindrical pore CD has radius $a$ and length $h$. The cylindrical surface CD of the pore has surface charge density $\sigma_{c}$, and the membrane surfaces BC and DE have surface charge density $\sigma_{m}$. An electrical potential difference is applied betwen the fluid reservoirs at either side of the membrane, and electroosmotic flow is generated by the resulting electric field acting on the charge cloud adjacent to the charged surfaces. The analysis of Mao et al.1 assumed that the external reservoirs on either side of the pore were unbounded, with radius $b=\infty$. For numerical work, the external reservoirs were bounded by uncharged cylinders of radius $b\gg a$, sufficiently large that numerical results when $h=0$ differed little from the analytic results for $h=0$ and $b$ infinite. There have been many studies in which flow is generated in cylinders of different dimensions, connected either in series 15 or in networks intended to represent porous media 16. Here, however, we are interested in the effect of the surfaces BC and DE of the membrane on electroosmotic flow within the cylindrical pore, and any boundaries AB, EF of the external reservoirs are so far away that they can be neglected.

In section 2 we set up the approximate analysis of end effects, and compare results to those obtained from full numerical computations. The analysis is presented from first principles, but can alternatively be set within the framework of the reciprocal theorem, as explained in section 2.6. The agreement between the approximate analysis and full computation is in general good, except for large Debye lengths $\kappa^{-1}\gg a$. In section 3 we consider this case in more detail, in order to evaluate how much of the charge cloud due to the charged walls of the cylindrical pore lies within the pore and how much spills out beyond the ends of the pore. When this overspill is taken into account, the agreement between the computations and the approximate model is improved.

2 2. Composite electroosmotic coefficient

2.1 2.1 The pore geometry

The axisymmetric geometry that we are considering is shown in 1. We use cylindrical polar coordinates $(r,z)$, with the $z$ axis along the axis of symmetry and $z=0$ at the midpoint of the cylindrical pore, the ends of which are at $z=\pm h/2$. When $h=0$ we shall also use oblate spherical coordinates $(\xi,\eta)$, with

[TABLE]

where $-\infty<\xi<\infty$ and $0\leq\eta<\pi/2$.

The cylindrical pore and the reservoirs at either end are filled with liquid with electrical conductivity $\Sigma$ and viscosity $\mu$. The wall CD of the cylindrical pore is charged, with uniform surface charge density $\sigma_{c}$, and the surface charge density over the membrane surfaces BC, DE, is $\sigma_{m}$. We assume that the reservoir boundaries AB, EF are uncharged and at infinity. We shall occasionally refer to the surface potential $\zeta$, which will not in general be uniform, but which is required to be small, with $\zeta\ll kT/e$, where $e$ is the elementary charge and $kT$ the Boltzmann temperature. The electrical potential $\phi_{0}$ within the equilibrium charge cloud therefore satisfies the linearized Poisson-Boltzmann equation, so that

[TABLE]

and the charge density in the equilibrium charge cloud is

[TABLE]

2.2 2.2 The applied electric field

The applied electric field is $\mathbf{E}=-\nabla\chi$, where the potential $\chi$ satisfies the Laplace equation

[TABLE]

with gradient

[TABLE]

normal to the walls of the membrane and of the cylindrical pore. In $z>0$, the electric potential far from the membrane is $\chi=\phi_{1}/2$, and the potential far from the membrane in $z<0$ is $\chi=-\phi_{1}/2$.

When the membrane thickness $h=0$, the potential can be expressed explicitly as 21

[TABLE]

On the plane of the membrane, within the circular opening,

[TABLE]

The liquid within the pore has electrical conductivity $\Sigma$; we have assumed that surface charge density (and hence the density of charge in the cloud of counter ions) is small, so that surface conductivity may be neglected. The total electric current $I_{m}$ flowing through the hole in the membrane is therefore

[TABLE]

If $h>0$, we assume that the potential within the cylindrical pore varies linearly and approximate the potential within the pore as

[TABLE]

as would be expected in the absence of any end effects. The potential in $z>h/2$ is approximated by that outside a membrane (with a hole) of zero thickness:

[TABLE]

with $\chi(r,z)=-\chi(r,-z)$. This approximation (9) and (10) is continuous at $z=\pm h/2$ where the potential is assumed to be $\phi_{2}$/2 across the entire width of the opening (by (7)). The as yet unspecified potential $\phi_{2}$ is determined by requiring continuity of the electrical current at $z=\pm h/2$. The current $I_{c}$ through the cylindrical pore is

[TABLE]

and the electrical current through the reservoir in $z>h/2$ is, by (8),

[TABLE]

Equating $I_{c}$ (11) and $I_{m}$ (12), we find

[TABLE]

2.3 2.3 Electroosmosis through an infinite cylindrical pore

We assume throughout this paper that the perturbation of the equilibrium charge cloud by the applied electric field and by fluid motion is negligibly small. The force acting on the ions in the charge cloud due to the applied electric field $-\nabla\chi$ is therefore $-\rho_{0}\nabla\chi$.

The equilibrium potential within an infinite cylindrical pore is

[TABLE]

In the absence of any end effects, if the electric field $E_{0}=-\phi_{2}/h$ is applied along the length of the cylindrical pore, the fluid velocity is 22

[TABLE]

and the total electroosmotic volumetric flow rate is 2

[TABLE]

where the electroosmotic coefficient

[TABLE]

2.4 2.4 Electroosmosis through a membrane ($h=0$)

It was shown by Mao et al. 1 that if the equilibrium charge density is $\rho_{0}$, the imposed electric field is $\mathbf{E}=-\nabla\chi$ and the fluid velocity generated by a pressure difference $p_{1}$ across a pore (of arbitrary geometry) is

[TABLE]

then the reciprocal theorem 23 for Stokes flows can be used to show that electroosmotically generated volumetric flow rate through the pore is

[TABLE]

where the integral is over all the fluid.

The fluid velocity generated by the pressure difference $p_{1}$ across a circular hole in a membrane of zero thickness is

[TABLE]

An explicit expression for $\mathbf{G}^{m}(r,z)$ is available 23, 1, and the potential $\chi$ is given by (6). The charge density in the equilibrium charge cloud around a membrane of zero thickness is 1

[TABLE]

which consists of the charge density adjacent to a uniform charged surface, from which has been subtracted the charge density around a uniformly charged disk. The integral (19) can be evaluated numerically 1, and the electroosmotic flow rate through a hole in a membrane of zero thickness can be expressed in the form

[TABLE]

where

[TABLE]

with

[TABLE]

2 shows a log-log plot of results for $H_{m}/H_{0}$ obtained by Mao et al. 1. The continuous line shows the analytic result (19) obtained via the reciprocal theorem, and the asymptote (23) for $a\kappa\ll 1$ is indicated.

The membrane has zero thickness, so that there is always a region near the edge of the pore where the Debye length $\kappa^{-1}$ cannot be considered small compared with $h$: Smoluchowski’s analysis for thin charge clouds, which would predict $H=6H_{0}/(a\kappa)$ if $\zeta_{m}$ took the uniform value $\epsilon\kappa\sigma_{m}$, therefore cannot automatically be invoked when $a\kappa\gg 1$. However, if we set up a local coordinate $s$ indicating distance from the edge of the pore, both the electric potential $\chi$ (6) and the fluid velocity $\mathbf{G}^{m}$ (20) vary as $s^{1/2}$ when $s\ll a$ (i.e. near the pore edge). The charge cloud density $\rho_{0}$ decays over a lengthscale $\kappa^{-1}$, and only counter-ions of membrane surface charge within a distance $\kappa^{-1}$ from the edge contribute to $\rho_{0}$ within the hole. The contribution of the edge to the integral (19) is therefore $O((a\kappa)^{-1})$, as was similarly found for the electrophoretic velocity of a charged disk 24. We therefore expect $H_{m}\sim KH_{0}/(a\kappa)$ when $a\kappa\gg 1$. The data in 2 do not extend to sufficiently high values of $a\kappa$ to allow us to estimate $K$ with any accuracy, and for the figure we simply indicate the line $K=6$ suggested by the Smoluchowski analysis. Clearly the asymptote corresponds to a value $K<6$. A similar reduction in the broadside electrophoretic velocity of a disk below the value predicted by Smoluchowski was noted by Sherwood & Stone.24 Individual points in 2 indicate results obtained from full numerical solutions of the Poisson-Nernst-Planck equations in a symmetric electrolyte at low applied potential and low surface charge. In the computations, the length of the reservoirs in the $z$ direction was equal to their radius $b$, with $b=\max(10a,10\kappa^{-1})$. Other details of the computations are reported in section 4.

2.5 2.5 Composite electroosmotic coefficient $H_{\rm comp}$

When $h>0$ it is natural to suppose that the electric field ouside the membrane pumps fluid towards the cylindrical pore at a rate

[TABLE]

and the electric field within the cylindrical pore pumps fluid through the pore at a rate

[TABLE]

However, in general, $Q_{me}$ (25) and $Q_{ce}$ (26) differ, and a pressure $\pm p_{2}/2$ builds up at $z=\pm h/2$ (i.e. at the entrance and exit to the cylindrical pore) in order to ensure that the volumetric flow rate is continuous. We now determine this pressure $p_{2}$.

Consider a membrane of zero thickness ($h=0$), with pressure $p=p_{1}/2$ (above the reference ambient pressure) at infinity on the side $z>0$, and with $p=-p_{1}/2$ at infinity on the other side. The pressure within the hole in the membrane is

[TABLE]

The fluid velocity generated by the pressure difference $p_{1}$ across the membrane is $\mathbf{u}=p_{1}\mathbf{G}^{m}$ (20), and the corresponding volumetric flow rate is 23

[TABLE]

If $h>0$ we approximate the pressure field in the fluid in much the same way as we approximated the electrical potential within the fluid: we patch a linearly varying pressure $p(z)$ within the cylindrical pore to the pressure field outside a membrane of zero thickness, and we take the pressure over the two ends $z=\pm h/2$ of the cylindrical pore to be $\pm p_{2}/2$. Thus the pressure within the pore is approximated as

[TABLE]

the fluid velocity within the pore is

[TABLE]

and the volumetric flow rate within the pore is

[TABLE]

Outside the cylindrical pore, the fluid velocity is assumed now to be

[TABLE]

with $u_{r}(r,z)=-u_{r}(r,-z)$ and $u_{z}(r,z)=u_{z}(r,-z)$. The volumetric flow rate outside the membrane is now

[TABLE]

We have ensured that the pressure (but not the fluid velocity nor the volumetric flow rate) is continuous across the ends $z=\pm h/2$ of the cylindrical pore.

When an electric field generates an electroosmotic velocity, the volumetric flow rates within the cylindrical pore and outside the membrane are identical if $p_{2}$ is such that $Q_{mh}+Q_{me}=Q_{ch}+Q_{ce}$, i.e. if

[TABLE]

But the pressure at infinity is zero in the electroosmotic problem, so $p_{1}=0$, and $\phi_{2}$ is given by (13). Hence

[TABLE]

and the total electro-osmotic flow is

[TABLE]

An alternative derivation of this approximate composite $H_{\text{comp}}$ (36) is given in the next section.

Inserting into (36) the various estimates for $G_{m}$ (33), $G_{c}$ (31), $R_{m}$ (8) and $R_{c}$ (11), we obtain

[TABLE]

For small $h/a$ the approximate composite $H_{\text{comp}}$ is larger than $H_{m}$ if

[TABLE]

Experimental arrangements often involve measurements at fixed current. Thus, it is helpful to define a coefficient that gives the electroosmotic flux per unit current. This quantity may be obtained readily from (11), (13) and (36):

[TABLE]

2.6 2.6 Composite electroosmotic coefficient $H_{\rm comp}$

derived via the reciprocal theorem

We now show that approximations to the electric potential $\chi$ and pressure-driven velocity $\mathbf{G}$ within a pore of non-zero length $h>0$, when inserted into the integral expression (19) for the electroosmotic volume flux, lead to an approximate electroosmotic coefficient identical to $H_{\text{comp}}$ (37) obtained in the previous section.

We have already shown that we may approximate the electric potential by a composite potential (9), (10), of the form

[TABLE]

We now create a similar approximation for the fluid velocity for flow through a membrane of thickness $h$ subjected only to a pressure drop $p_{1}$ but no applied potential drop. We suppose that in $z>h/2$ the fluid velocity is given by (32), corresponding to flow outside a membrane of zero thickness, and that within the cylindrical pore the fluid velocity is given by (30). Continuity of the volumetric flow rates (31), (33) at the entrance to the cylindrical pore requires that the pressure $\pm p_{2}/2$ at the two ends of the pore satisfies

[TABLE]

so that

[TABLE]

Hence our approximation to the fluid velocity is $\mathbf{u}=\mathbf{G}p_{1}$, with

[TABLE]

We now use the approximations (40) and (43) in the integral (19) in order to compute the electroosmotic volumetric flow rate. But the integration splits naturally into an integral over the cylindrical pore and an integral over the regions outside the membrane. The integral over the cylindrical pore is exactly the integral required to determine the electroosmotic flow rate $H_{c}$ (16) in a cylinder, and the integral outside the membrane is exactly that required to determine $H_{m}$ (22). Hence the integral yields the composite electroosmotic flow rate

[TABLE]

identical to (36), obtained in section 2.5 by elementary methods.

2.7 2.7 Predictions of the composite electroosmotic coefficient

3 shows $H_{\text{comp}}$ (37) as a function of $h/a$, for four different values of $a\kappa$, with $\sigma_{m}=\sigma_{c}$. Also shown are the results of full numerical computations based on the Poisson-Nernst-Planck equations 1. We see that for $a\kappa\geq 1$ the approximate analysis captures the main features of the full numerical results, and it is clear from (37) that it also has the correct limits as $h/a\rightarrow 0$ and $h/a\rightarrow\infty$. However, it is also evident from 3(d) that the theory is unsatisfactory when $a\kappa\ll 1$. We discuss this limit in the next section, where we shall show that when $a\kappa\ll 1$ some of the charge cloud of ions that neutralizes the surface charge on the cylindrical wall of the pore spills out of the ends of the pore, where it is less effective at generating electroosmotic flow. The scenario is shown schmatically in 4.

3 3. Charge overspill from the ends of the pore, $a\kappa\ll 1$

3.1 3.1. Overspill of charge from the end of a semi-infinite pore

We consider a cylindrical pore of radius $a$, with surface charge density $\sigma_{c}$. When the Debye length $\kappa^{-1}\gg a$, the equilibrium potential $\phi_{0}$ (14) in an infinitely long cylinder can be expanded as

[TABLE]

where

[TABLE]

Thus the equilibrium potential $\phi_{0}$ and charge density $\rho_{0}=-\epsilon\kappa^{2}\phi_{0}$ within the charge cloud vary little over the cross-section of the pore. On the other hand, if the cylinder is not infinitely long and uniform, $\phi_{0}$ and $\rho_{0}$ vary in the axial ($z$) direction with a length scale $\kappa^{-1}$. We can therefore consider the equilibrium potential $\phi_{0}$ within the cylindrical pore to be a function only of $z$ 25.

We first consider a semi-infinite, charged cylindrical pore going from $z=0$ to $z=\infty$. The equilibrium potential $\phi_{0}$ satisfies a one-dimensional Poisson-Boltzmann equation,

[TABLE]

The solution that tends to the uniform potential $\phi_{a}$ within the pore as $z\rightarrow\infty$ far from the pore end at $z=0$, is

[TABLE]

for some unknown constant $A$. The charge density within the charge cloud inside the pore is $-\epsilon\kappa^{2}\phi_{0}$, and when the cylindrical pore is infinite (and hence uniform) the charge per unit length in the charge cloud is $-\pi a^{2}\epsilon\kappa^{2}\phi_{a}=-2\pi\sigma_{c}a$, equal and opposite to the charge per unit length on the pore walls. When the pore is semi-infinite, with a non-uniform charge cloud (48), the total charge that is lost from within the pore is

[TABLE]

At the end of the pore ($z=0$) the potential is $\phi=\phi_{a}-A$.

In $z<0$ the charge cloud is no longer confined by the walls of the cylindrical pore and spreads out radially: it is no longer possible to assume that $\phi_{0}$ is a function of $z$ alone. We therefore need to solve the linearized Poisson-Boltzmann equation in the half-space $z<0$, with $\phi_{0}=\phi_{a}-A$ over the region $z=0$, $r<a$ and $\partial\phi_{0}/\partial z=0$ on $z=0$, $r>a$. At large distances from the end of the pore, the potential decays as $\exp(-\kappa R)/R$, where $R=(z^{2}+r^{2})^{1/2}$ is a spherical polar coordinate, but in the important region $R=O(a)$ the potential can be approximated by the electrostatic potential corresponding to a solution of the Laplace equation (i.e. $\kappa=0$). Hence, from (6),

[TABLE]

To relate the potential (50) to the amount of charge in the overspilling charge cloud (in $z<0$), we note that the charge on one side of a charged disk at uniform potential $(\phi_{a}-A)$ in unbounded space is $q=4a\epsilon(\phi_{a}-A)$. Alternatively, one can argue that far from the plane $z=0$, the spherical distance $R\approx a\cosh\xi$, so that the potential (50) is approximately

[TABLE]

In a spherically symmetric geometry this field corresponds to the far field around a point charge of magnitude $8a\epsilon(\phi_{a}-A)$, and the total surface charge on one side of the disk is $q=4a\epsilon(\phi_{a}-A)$, in agreement with the charge obtained by considering the capacitance of the disk. The charge in the overspill charge cloud in $z<0$ is equal and opposite to $q$, and is therefore

[TABLE]

But the total charge in the overspill outside the end of the pore must be equal and opposite to the charge that has been lost from within the pore. Hence

[TABLE]

so that

[TABLE]

and the potential at the end of the pore is

[TABLE]

The charge that has been lost from the end of the pore is equivalent to the charge usually found in a pore of length

[TABLE]

3.2 3.2. Overspill from the two ends of a finite pore

We can now perform the same analysis for a pore that occupies the region $-h/2<z<h/2$. The equilibrium potential within the pore has the form

[TABLE]

where we have chosen the solution that is symmetric about the centre of the pore at $z=0$. The charge that has been lost from within the pore is

[TABLE]

The total flux of electric field through the two ends of the pore is

[TABLE]

Comparing (58) and (59) we conclude that $C=0$. The potential over the ends of the pore is

[TABLE]

The total charge in the two overspill charge clouds is therefore, by (52),

[TABLE]

and this must be equal to the charge (58) lost from within the pore. Hence

[TABLE]

and

[TABLE]

The total charge that has been lost (from the two ends) is equivalent to a total lost length

[TABLE]

We see from eqs (56) and (65) that when $\kappa h\gg 1$ the lost charge is twice that lost from a single end of a pore. We also note that $h-h_{\text{lost}}>0$, and that

[TABLE]

3.3 3.3. Overspill from the membrane surface into the pore

If the cylindrical pore itself is uncharged, but the membrane surfaces are charged, ions from the charge cloud adjacent to the membrane surface are able to move into the ends of the pore.

If the membrane has zero thickness, the charge density $\rho_{0}$ in the equilibrium charge cloud is given by (21), and both $\rho_{0}$ and the potential $\phi_{0}=-\rho_{0}/(\epsilon\kappa^{2})$ vary over the area of the pore. Nevertheless, we may work out the mean potential over the circular pore,

[TABLE]

where, when $a\kappa\ll 1$,

[TABLE]

Thus when the membrane has zero thickness (and there is no cylindrical pore into which ions can escape) the absence of surface charge over the area of the pore changes the average potential over the opening from the value $\phi_{0}=\sigma_{m}/(\epsilon\kappa)$ due to a uniformly charged surface to $\beta\sigma_{m}/(\epsilon\kappa)$, where

[TABLE]

We now consider the charge that leaks into a pore of length $h>0$ from the charge clouds on either side of the membrane. We suppose that the potential on the planes $z=\pm h/2$ is perturbed by an amount $D$, and becomes

[TABLE]

Within the pore, the potential obeys the 1-dimensional Poisson-Boltzmann equation (47), with solution

[TABLE]

and the additional charge within the pore is

[TABLE]

Outside the pore, the perturbed potential (71) is associated with a total additional charge (52)

[TABLE]

on the two sides of the membrane. But the total change in charge caused by this redistribution must be zero, i.e. $q_{\text{in}}+q_{\text{out}}=0$. Hence

[TABLE]

i.e.

[TABLE]

The total charge $q_{\text{in}}=-q_{\text{out}}$ (74) that leaks into the pore at the two ends corresponds to the charge inside a uniformly charged cylinder with surface charge density $\sigma_{m}$, of length

[TABLE]

Thus $h_{\text{gained}}$ (77) is smaller than $h_{\text{lost}}$ (64) by a factor $a\kappa\beta/2$. We can compare the predictions (64) and (77) against results obtained from full numerical solution of the nonlinear Poisson-Boltzmann equation with either $\sigma_{m}=0$ and $ae\sigma_{c}/(\epsilon kT)=a\kappa e\zeta_{c}/(kT)=0.00273$, or $\sigma_{c}=0$ and $ae\sigma_{m}/(\epsilon kT)=0.00273$: results for $a\kappa=0.1$ are given in Table 1. We see that there is excellent agreement between the numerical computations and the analysis presented above.

3.4 3.4. Composite electroosmotic coefficient

We first consider how the electroosmotic coefficients $H_{c}$ and $H_{m}$ are modified by the overspill of the charge cloud from inside the cylindrical pore to outside the membrane. If a uniform electric field of strength $E=-\phi_{1}/h$ is applied between the ends of the pore, the Navier Stokes equations for steady flow yield the axial velocity profile

[TABLE]

so that the volumetric flow rate is

[TABLE]

But $Q$ is independent of $z$ (by incompressibility), and the pressure $p$ difference between the two ends of the capillary is zero. Hence, integrating (81) along the length $h$ of the cylindrical pore, and noting that the total amount of charge in the charge cloud remaining within the pore is $2\pi a\sigma_{c}(h-h_{\text{lost}})$, we find

[TABLE]

which may be compared to the result (LABEL:H_c_kappa_small) which ignores overspill. The charge cloud outside the pore is enhanced by the overspill, and becomes (in $z>0$)

[TABLE]

with the final term $[\phi_{a}-B\cosh(\kappa h/2)]$, corresponding to the overspill charge cloud (63), being approximately valid in a volume $O(a^{3})$ around the pore, but invalid at large distance $O(\kappa^{-1})$ from the pore, where the exponential decay of the charge density is not captured by the solution (50) of the Laplace equation. The volumetric flow rate through a pore of zero thickness created by a potential difference $\phi_{1}$ is given by the integral (19) and was shown by Mao et al. 1 to be:

[TABLE]

where

[TABLE]

Hence the electroosmotic flow rate $Q=H_{m}\phi_{1}$ due to the charge cloud outside the membrane is modified, and $H_{m}$ becomes

[TABLE]

If $\sigma_{m}$ is comparable to $\sigma_{c}$, then we saw in Section 3.3 that the change in the charge within the pore due to the charge cloud outside the membrane entering the pore is $O(a\kappa)$ smaller than the loss of charge from the charge cloud within the pore to the regions outside the membrane. However, this contribution can be included with very little effort, and becomes important in the limit $h\rightarrow 0$, when the gain (78) in charge within the pore from the outside surface charge density $\sigma_{m}$ is proportional to $h_{\text{gained}}\propto h$, whereas the charge cloud (due to $\sigma_{c}$ within the pore) remaining within the pore, is proportional to $h-h_{\text{lost}}\propto h^{2}$, by (67). The electroosmotic coefficient $H_{c}$ for the cylindrical pore (82) becomes

[TABLE]

and the electroosmotic coefficient $H_{m}$ for the charge cloud outside the membrane (86) becomes

[TABLE]

Now that $H_{c}$ (87) and $H_{m}$ (88) have been corrected for the effects of overspill in the two directions, they can be inserted into the expression (37) for the composite electroosmotic coefficient $H_{\text{comp}}$. Results are shown in 5, together with full numerical solutions of the Poisson-Nernst-Planck equations. We see that the agreement between theory and computation is much better than when overspill is ignored (3(d)). Charge overspill or underspill causes the total charge of mobile ions within the pore to be different from what might be expected on the basis of net electroneutrality of the pore. Thus, the driving force is modified leading to deviations from the calculated result that ignores such effects. The “lost length” $h_{\text{lost}}$ in (87) restores this effect.

Note that when $h\ll\kappa^{-1}$ the effective length of the cylindrical pore $h-h_{\text{lost}}\approx\pi a\kappa^{2}h/2$, by (67). The approximation (87) for $H_{c}$ is therefore dominated by the term $h_{\text{gained}}$, and gives $H_{c}\sim\pi a^{4}\kappa\beta\sigma_{m}/(8h\mu)$, with $H_{c}/H_{m}\sim 3\pi a\beta/(8h)$. We conclude from (38) that $H_{\text{comp}}$ is a decreasing function of $h$ near $h=0$, as seen in 5.

4 4. Numerical simulation

We solve the time–independent PNP–Stokes equations numerically using the finite volume method. The governing equations are:

[TABLE]

Here $n^{i}$ is the local concentration of the i–th ion, $i=1\sim N$ and $N=2$ in this case, $z_{i}$ is the valence of the i–th ion, and in our case we have $z_{1}=-z_{2}=1$, $\omega^{i}$ is the mobility of the i–th species of ion. The electric potential is $\phi$, the flow velocity is $\mathbf{u}$ and $p$ is the pressure.

We solve the system of coupled equations (89)-(92) in the axisymmetric geometry depicted in 1. Thus, we have a cylindrical pore of radius $a$ and length $h$ connecting two large cylindrical reservoirs of radius $b$. The lengths of AB and EF in our simulation are also taken to be $b$. Numerically, $b$ is kept much larger than either $a$ or the Debye length $\kappa^{-1}$ so that the reservoirs are effectively infinite.

We apply the following boundary conditions. At A and F, the two ends of the reservoirs, ion concentrations are made equal to the concentration in the bulk electrolyte, or $n^{i}=n^{i}_{\infty}$; a potential difference of $\Delta V$ is applied across the system by setting $\phi$ to $\pm\Delta V/2$ respectively at A and F while the pressure at these two ends are set equal to the bulk pressure, $p=p_{\infty}$. At AB and EF, the side walls of the cylindrical reservoir, the radial component of the electric field, ionic flux and velocity are all set to zero, as the walls are considered to be far away from the pore. A zero tangential shear stress is imposed on AB and EF as well. At the membrane and pore surfaces BC, CD and DE, a no–flux condition is used for (90), and a no–slip condition is used for the flow. At solid fluid interfaces the electric potential is continuous but the normal component of the electric field undergoes a jump $[\epsilon\mathbf{E}\cdot\hat{n}]=\sigma_{j}$, where $j=m$ at BC and DE, $j=c$ at CD and we have used square brackets to indicate change of the enclosed quantity across the interface. The suffixes $m$ and $c$ indicate the surface charge densities on the face of the membrane and within the pore respectively. Here $\hat{n}$ is the unit normal at the surface directed into the fluid.

An electrohydrodynamic solver was implemented to solve the system described above using the OpenFOAM CFD library 26, a C++ library designed for computational mechanics. Structured mesh is constructed with polyMesh, a meshing tool built in OpenFOAM. The grid is refined near the membrane and pore surfaces to resolve the Debye layer. Grid independence was checked in all cases by refining the grid and verifying that the solution does not change within specified tolerances.

For the finite volume discretization of the governing equations, central differences are used for all diffusive terms in (90) and viscous terms in (91). A second–order upwind scheme is used for the convective terms in (90). The discretized linear system is solved using a pre-conditioned conjugate gradient solver if the matrix is symmetric or a pre-conditioned bi–conjugate gradient solver if the matrix is asymmetric 27.

An iterative scheme is used to solve the PNP–Stokes equations. Initially, the flow velocity is set to zero. Equations (89) and (90) are then solved sequentially in a loop with under-relaxation (to ensure stability of the nonlinear PNP system) until the absolute residual is smaller than a specified tolerance – in our case, $10^{-6}$. The electric force density $-\nabla\phi\sum_{i}z_{i}en^{i}$ is then obtained from this solution and used as an explicit external forcing in the next step - the solution of the incompressible Stokes flow problem: (91) and (92). The SIMPLE algorithm is used to solve the hydrodynamic problem. The flow field so computed is then substituted into (90) and the PNP equations are then solved again using the updated flow field. An outer loop is constructed to iterate over the PNP loop and Stokes flow module, until the solution changes negligibly between two outer iterations.

Our main object of interest is the volumetric flux, $Q$. This is obtained in a post processing step by numerically integrating the axial velocity over the plane $z=0$. At the low voltages employed, a substantially linear relation is found between $Q$ and $\Delta V$, the electroosmotic coefficient $H$ is then obtained from $H=Q/\Delta V$. The quantities $h_{lost}$ and $h_{gained}$ were obtained numerically by integrating the local charge density within the confines of the pore.

5 4. Concluding remarks

The analysis presented here shows that it is possible to use simple analyses based on continuity of volumetric flow rate and electric current to estimate electroosmotic end effects in a charged cylindrical pore traversing a membrane of thickness $h>0$. Note that we have made repeated use of the assumption that surface charge densities, and corresponding zeta potentials, are small. Not only have we worked with the linearized Poisson-Boltzmann equation (2), but we have used superposition to combine various contributions to the charge clouds due to overspill of the clouds from one region (inside/outside the pore) to the other. At high potentials it would also be necessary to keep track of the fluxes of individual ion species, rather than simply ensuring that the total electrical current is continuous 15.

The assumption of small potentials also justifies the neglect of other nonlinear electrokinetic effects such as Induced Charge Electroosmosis (ICEO) 17, 18 which could produce vortices in the vicinity of sharp corners 19 or near rapid constrictions in channels 20. However, numerical solutions confirm the expectation that the flow rate is only weakly affected by such vortices, particularly under conditions of small potentials 13

Recent experimental work 8, 9, 10, 11, 12 on nanopores report potential differences of $\Delta\phi\sim 0-200$ mV applied across the pore. Here we have assumed that $\Delta\phi\ll\zeta$, where $\zeta$ itself is assumed small in comparison with the thermal voltage $kT/e\sim 25$ mV. Thus, our results can only be expected to describe the initial linear part of the current-voltage and flow-voltage characteristics, even though numerical simulations seem to show 13 that this linear regime extends to applied voltages $\sim 100$mV.

Finally, we point out that the correction factor $\beta$ (70) reminds us that the hole in the charged membrane removes a circular region of surface charge and reduces the equilibrium potential at the entrance to the pore. The introduction of $\beta<1$ improved the agreement between theoretical and numerical results for $h_{\text{gained}}$ in Table 1. However, the analysis is not rigorous, since the equilibrium potential across the hole is not uniform. The $O(1-\beta)$ correction to the equilibrium potential corresponds to an $O(1-\beta)$ correction to the charge density $\rho_{0}$. If we use this in the integral expression (84) in order to determine a correction to the electroosmotic flow rate through a membrane of zero thickness, the analysis suggests that the correction to the leading order result (23) for $a\kappa\ll 1$ should be $O((a\kappa)^{2})$, whereas investigation of the difference (seen in 2) between numerical results and the asymptote (23) indicates additional corrections $O((a\kappa)^{2}\ln a\kappa)$.

{acknowledgement}

JDS thanks the Department of Applied Mathematics, University of Cambridge, and the Institut de Mécanique des Fluides de Toulouse, for hospitality. MM & SG acknowledge support from the NIH through Grant 4R01HG004842.