# Analytic theory of nonlinearly coupled electrokinetics in nanochannels

Yuki Uematsu (植松 祐輝)1,2<EMAIL_ADDRESS>1Department of Physics and
Information Technology, Kyushu Institute of Technology, 820-8502 Iizuka,
Fukuoka, Japan
2PRESTO, Japan Science and Technology Agency, Kawaguchi, Saitama 332-0012,
Japan

###### Abstract

The nonlinear electrokinetic response of ionic solutions is important in
nanofluidics. However, quantitatively understanding the mechanisms is still a
challenging problem because of a lack of analytic approaches. Here, a general
framework for calculating the nonlinear electrokinetic coefficients of
strongly confined electrolytes is constructed using a perturbation scheme of
the pressure and voltage differences across a nanochannel. The theory is
applied to an electrically neutral nanochannel filled with electrolytes, and
analytic expressions for the first- and third-order electrokinetic
coefficients are obtained. We demonstrate that the combination of high
hydrodynamic permeability and ion-wall friction plays an essential role in
nonlinear electrokinetics. Furthermore, we analytically demonstrate that the
external flow induces uniform excess charge inside the nanochannel.

##  I Introduction

Recent developments in the technology to manipulate, manufacture, and
fabricate atoms and molecules at the nanoscale level enable the production of
nanofluidic devices that exhibit excellent transport and electrokinetic
properties. [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16] One of the
examples is the ultra-high hydrodynamic permeability of water flowing through
carbon nanotubes. [3, 4, 5] The carbon surface was experimentally demonstrated
to be hydrophobic, with the slip length on the micrometer scale. [3, 4, 5] A
large zeta potential [7] and anomalous conductance [8] of electrolytes
confined in nanochannels are also useful for achieving energy conversion at
high efficiency. In addition to these linear responses, electrolyte solutions
in strong confinements, typically sub- or few nanometers, exhibit unique
nonlinear electrokinetics. [7, 9, 10, 15] In strong confinements, the
streaming current depends on the applied voltage and the applied pressure
difference. More precisely, the ionic current response has a component
proportional to $(\Delta p)(\Delta\psi)^{2}$ and $(\Delta p)^{2}(\Delta\psi)$,
where $\Delta p$ is the applied pressure difference, and $\Delta\psi$ is the
applied voltage difference [7, 10]. This unique non-linearity is essentially
different from the previously reported nonlinear electrokinetics [17, 18, 19]
because in previous studies, the focus was on the electro-osmotic/phoretic
velocity proportional to $(\Delta\psi)^{2}$ and $(\Delta\psi)^{3}$, where
$\Delta p$ and $\Delta\psi$ are not nonlinearly coupled. Janus particles
exhibit a quadratic dependence of the electrokinetic flow
$\sim(\Delta\psi)^{2}$ because the induced charge on the particle surface is
linear with $\Delta\psi$. [20] The cubic dependence of the electrophoretic
mobility is a result of strong surface conduction in the spherical geometry.
[21, 18]

Phenomenological nonlinear polynomials for electrokinetic flux were first
suggested for glass membranes, ion-exchange membranes, and animal membranes.
[22] However, these experimental data exhibit complex dependence rather than
simple polynomial nonlinear dependence of the applied pressure and voltage.
Nevertheless, the recent nanofluidic measurements on the streaming current [7,
10] excellently fit parabolic curves represented by $\sim(\Delta p)^{2}$ or
$\sim(\Delta\psi)^{2}$, implying that there should exist a novel fundamental
physical mechanism of the nonlinear responses. At present, numerical
simulation of electro-hydrodynamic equations including high hydrodynamic
permeability, ion-wall friction, and energy barriers inside the nanochannels
can reproduce the quadratic responses of the streaming current on the applied
pressure or the voltage, [7, 10] but the analytic calculation of each
nonlinear electrokinetic coefficient is still lacking.

Here, we construct a general framework to solve electrohydrodynamic equation
by perturbation theory with the applied pressure and voltage difference. By
applying this framework to a confined electrolyte solution in an uncharged
nanochannel, the linear and nonlinear electrokinetic coefficients are
analytically calculated. The obtained analytic expression has a simple form in
the limit of a long channel, and these equations demonstrate that the
combination of the high hydrodynamic permeability and ion-wall friction is the
necessary condition for experimentally measurable nonlinear responses.
Furthermore, we investigate how the electrostatic potential profiles along the
channel coordinate respond to the perturbations in the applied voltage and
water flow. Calculation of the induced charge density demonstrates that water
flow induces uniform effective charges that are quadratic with the applied
voltage and the water flow. The flow-induced charge inside the nanochannel is
a crucial factor of the third-order electrokinetics in our model. In the end,
we discuss the application of the theory to the previous experiments. [10]
Although all of the quantities in the theory cannot be determined by the
limited experimental data, the single carbon nanotube exhibits the
dimensionless nonlinear electrokinetic coefficients in the order of $10^{3}$
to $10^{5}$, which is significantly large. Such large nonlinear
electrokinetics suggests the combination of large hydrodynamic permeability
and ion-wall friction in the nanochannels. Furthermore, the joint length
between the nanochannel and reservoirs is found to be important to quantifying
the nonlinear electrokinetic coefficients, whereas the experimental definition
of the joint length has not been investigated in the previous nanochannel
experiments.

Figure 1:  Nanochannel connected with reservoirs filled with monovalent
electrolyte solution. The pressure difference, $\Delta p$, and the voltage
difference, $\Delta\psi$, can be externally applied, whereas the solution
flux, $Q$, and the ionic current, $I$, are measurable.

##  II Model

###  II.1 1-dimensional electrohydrodynamic equations

The model is similar to those suggested to consider the confinement effect in
nanochannels [7, 10]. We consider an uncharged nanochannel connected with two
reservoirs filled with electrolyte solutions as illustrated in Fig. 1. The
system is assumed to be axisymmetric with respect to the $x$-axis, and the
surface of the nanochannel and reservoirs is given by $R(x)$. In the
nanochannel region, $R(x)=R_{0}$ is constant, where $R_{0}$ is the radius of
the nanochannel, and $L$ is the length of the nanochannel. The electrolyte is
monovalent and fully dissociated in the solution with the dielectric constant
$\varepsilon$ and the viscosity $\eta$. The electrolyte concentration in the
reservoirs is set at $\rho_{0}$, and the system is stored at the temperature
$T$. The cations and anions are assumed to have the same diffusion constant
$D$ in the reservoirs and the nanochannel, but we can easily release this
assumption. The inhomogeneity of the ion concentration, the electrostatic
potential, and ion and water velocity fields in the radial and angular
coordinates are neglected, and their inhomogeneity is only considered in the
axial direction. Furthermore, the ionic volume interactions [23] are neglected
even though the salt concentration is high, and anisotropy and reduction of
the dielectric constant inside the nanochannel [24] are neglected. In the
nanotube, the ions feel the energy barrier due to the different solvation
structure from bulk, but this is omitted because we consider the high salinity
where the electrostatic screening reduces the energy barrier [25]. These
simplifications yield a set of electro-hydrodynamic equations for the pressure
field $p(x)$, the water velocity field $v_{\mathrm{w}}(x)$, the ion velocity
field $v_{i}(x)$ ($i=\pm$ denotes cation and anion), the $i$ type ion
concentration $\rho_{i}(x)$, and the electrostatic potential $\psi(x)$ given
as follows:

$\displaystyle-\frac{d}{dx}\left[\pi R(x)^{2}v_{\mathrm{w}}(x)\right]=0,$ (1)
$\displaystyle-\frac{d}{dx}\left[\pi R(x)^{2}\rho_{i}(x)v_{i}(x)\right]=0,$
(2)
$\displaystyle-\frac{dp(x)}{dx}-\sum_{i=\pm}\rho_{i}(x)\xi_{i}(x)[v_{\mathrm{w}}(x)-v_{i}(x)]\qquad$
$\displaystyle-\frac{2\lambda_{0}(x)}{R(x)}v_{\mathrm{w}}(x)=0,$ (3)
$\displaystyle-\frac{d\mu_{i}}{dx}-\xi_{i}(x)[v_{i}(x)-v_{\mathrm{w}}(x)]-\lambda_{i}(x)v_{i}(x)=0,$
(4)
$\displaystyle\frac{1}{R(x)^{2}}\frac{d}{dx}\left[\varepsilon\varepsilon_{0}R(x)^{2}\frac{d}{dx}\psi(x)\right]=-e\sum_{i=\pm}q_{i}\rho_{i}(x),$
(5)

where $\xi_{i}(x)$ is the friction coefficient between ion and water,
$\lambda_{0}(x)$ is the friction coefficient between water and wall per
surface area, $\lambda_{i}(x)$ is the friction coefficient between ion and
wall, $\varepsilon_{0}$ is the vacuum permittivity, $e$ is the elementary
charge, and $q_{i}(=\pm 1)$ is the valence of the ions. The electrochemical
potential of ions, $\mu_{i}$ is defined by

$\mu_{i}(x)=k_{\mathrm{B}}T\ln\frac{\rho_{i}(x)}{\rho_{0}}+eq_{i}\psi(x),$ (6)

where $k_{\mathrm{B}}$ is the Boltzmann constant. Eqs. 1 and 2 are the
continuity equations for the water flux and the ionic flux. Eqs. 3 and 4 are
the force balance equations for water and ions, whereas eq. 5 is the Poisson
equation. A similar model was first introduced for explaining nonlinear
electrokinetics in the subnanometer channels [7] and also in a single carbon
nanotube. [10] Eqs. 3 and 4 can be derived by the Onsager variational
principle [26] as described in supplementary material.

For analytical tractability, we introduce the profile of the ionic diffusion
constant $D_{i}(x)$ and the profile of the ionic friction ratio $\beta_{i}(x)$
as follows:

$\displaystyle D_{i}(x)$ $\displaystyle=$
$\displaystyle\frac{k_{\mathrm{B}}T}{\xi_{i}(x)+\lambda_{i}(x)},$ (7)
$\displaystyle\beta_{i}(x)$ $\displaystyle=$
$\displaystyle\frac{\xi_{i}(x)}{\xi_{i}(x)+\lambda_{i}(x)}.$ (8)

Furthermore, in this study, we consider a specific geometry $R(x)$ explicitly
given by:

$R(x)=\left\\{\begin{array}[]{lcc}\displaystyle{R_{0}}&\mathrm{for}&|x|<\displaystyle\frac{L}{2}\\\
\\\
\displaystyle{R_{0}}\mathrm{e}^{\frac{|x|-L/2}{2\delta}}&\mathrm{for}&|x|>\displaystyle\frac{L}{2}\\\
\end{array},\right.$ (9)

where $R_{0}$ is the radius of the nanochannel, $L$ is the nanochannel length
and $\delta$ is the characteristic length at the joint between the nanochannel
and reservoirs. In addition, we use a sharp step function for $\lambda_{0}(x)$
and $\beta_{i}(x)$ given by:

$\displaystyle\lambda_{0}(x)$ $\displaystyle=$
$\displaystyle\frac{\eta}{b}\theta(\frac{L}{2}-|x|),$ (10)
$\displaystyle\beta_{i}(x)$ $\displaystyle=$
$\displaystyle\alpha_{i}+(1-\alpha_{i})\theta(|x|-\frac{L}{2}),$ (11)

where $\eta$ is the solution viscosity, $b$ is the slip length of the
nanochannel, and $\alpha_{i}$ is the ion-wall friction ratio. The case of
$\alpha_{i}=1$ corresponds to no ion-wall friction force, and the case of
$\alpha_{i}=0$ corresponds to only ion-wall friction (no ion-water friction).
Furthermore, we use the homogeneous diffusion profile $D_{+}(x)=D_{-}(x)=D$
for both cations and anions. This limitation can be easily released in the
following analytical calculation. The boundary conditions of the equations
are:

$\displaystyle p(-H)$ $\displaystyle=$ $\displaystyle-p(H)=\frac{\Delta
p}{2},$ (12) $\displaystyle\psi(-H)$ $\displaystyle=$
$\displaystyle-\psi(H)=\frac{\Delta\psi}{2},$ (13) $\displaystyle\rho_{i}(-H)$
$\displaystyle=$ $\displaystyle\rho_{i}(H)=\rho_{0},$ (14)

where $H$ is half of the system size, $\Delta p$ is the applied pressure
difference, and $\Delta\psi$ is the applied voltage.

###  II.2 Derivation of the integro-differential Poisson equation

In this section, we derive an integro-differential Poisson equation, eq. 27
from eqs. 1 to 5. The similar method was first introduced in Ref. 10. Eqs. 1
and 2 are directly solved by introducing the three integral constants: the
water flux, $Q$, and the $i$ type ionic current, $J_{i}$. The water velocity
and ionic currents are solved as:

$\displaystyle v_{\mathrm{w}}(x)$ $\displaystyle=$ $\displaystyle\frac{Q}{\pi
R^{2}(x)},$ (15) $\displaystyle\rho_{i}(x)v_{i}(x)$ $\displaystyle=$
$\displaystyle\frac{J_{i}}{\pi R^{2}(x)}.$ (16)

At first, we introduce the effective potential, $W_{i}(x)$, and the effective
concentration, $c_{i}(x)$, as

$\displaystyle W_{i}(x)$ $\displaystyle=$ $\displaystyle
q_{i}\varphi(x)-\mathrm{Pe}\frac{\pi
R_{0}^{2}}{L}\int^{x}_{0}\frac{\beta_{i}(x^{\prime})}{\pi
R(x^{\prime})^{2}}dx^{\prime},$ (17) $\displaystyle c_{i}(x)$ $\displaystyle=$
$\displaystyle\frac{\rho_{i}(x)}{\rho_{0}}\mathrm{e}^{W_{i}(x)},$ (18)

where the Péclet number $\mathrm{Pe}$ is defined by:

$\mathrm{Pe}=\frac{QL}{\pi R_{0}^{2}D}.$ (19)

Solving eq. 4 with respect to the ionic velocity $v_{i}$, and plugging eqs. 17
and 18 into $v_{i}$, we obtain:

$v_{i}(x)=-\frac{1}{\xi_{i}(x)+\lambda_{i}(x)}\frac{d\mu_{i}}{dx}+\beta_{i}(x)v_{\mathrm{w}}(x)=-\frac{D}{c_{i}(x)}\frac{dc_{i}}{dx}.$
(20)

Because $J_{i}=\pi R^{2}(x)\rho_{i}(x)v_{i}(x)=\textrm{const.}$, we obtain:

$\frac{dc_{i}}{dx}=-\frac{J_{i}}{\rho_{0}D}\frac{\mathrm{e}^{W_{i}(x)}}{\pi
R(x)^{2}}.$ (21)

Eq. 21 can be solved with the boundary conditions $\rho_{i}(-H)=\rho_{0}$ and
$\varphi(-H)=\Delta\varphi/2$, as:

$c_{i}(x)=\mathrm{e}^{\frac{q_{i}\Delta\varphi+p_{i}\mathrm{Pe}}{2}}-\frac{J_{i}}{\rho_{0}D}f_{i}(x),$
(22)

where

$\displaystyle p_{i}$ $\displaystyle=$ $\displaystyle\frac{\pi
R_{0}^{2}}{L}\int^{H}_{-H}\frac{\beta_{i}(x^{\prime})dx^{\prime}}{\pi
R(x^{\prime})^{2}},$ (23) $\displaystyle f_{i}(x)$ $\displaystyle=$
$\displaystyle\int^{x}_{-H}\frac{\mathrm{e}^{W_{i}(x^{\prime})}}{\pi
R(x^{\prime})^{2}}dx^{\prime}.$ (24)

The boundary conditions, $\rho_{i}(H)=\rho_{0}$ and
$\varphi(H)=-\Delta\varphi/2$, yield

$c(H)=\mathrm{e}^{W_{i}(H)}=\mathrm{e}^{-(q_{i}\Delta\varphi+p_{i}\mathrm{Pe})/2},$
(25)

and solving eqs. 22 and 25 with respect to $J_{i}$ gives

$J_{i}=\frac{2\rho_{0}D}{f_{i}(H)}\sinh\frac{q_{i}\Delta\varphi+p_{i}\mathrm{Pe}}{2}.$
(26)

Therefore, eqs. 1 to 5 can be converted into a single integro-differential
Poisson equation rewritten as:

$\frac{1}{\kappa^{2}R(x)^{2}}\frac{d}{dx}\left[R(x)^{2}\frac{d}{dx}\right]\varphi=-\frac{1}{2}\sum_{i=\pm}q_{i}c_{i}(x)\mathrm{e}^{-W_{i}(x)},$
(27)

where $\kappa^{2}=2e^{2}\rho_{0}/\varepsilon\varepsilon_{0}k_{\mathrm{B}}T$,
$c_{i}(x)$, $f_{i}(x)$, and $W_{i}(x)$ are defined by eqs. 17, 22, 24, and 26.

For given $Q$ (or $\mathrm{Pe}$) and $\Delta\varphi$, eq. 27 yields the
electrostatic potential profile $\varphi(x)$. The individual ionic current,
$J_{i}$, can be calculated by eqs. 24 and 26, and the ionic current is given
by:

$I=e\sum_{i}q_{i}J_{i}.$ (28)

The pressure distribution can be calculated by integrating eq. 3 as:

$\begin{split}p(x)=&-\int^{x}_{0}\frac{2\lambda_{0}(x^{\prime})}{R(x^{\prime})}\frac{Q}{\pi
R(x^{\prime})^{2}}dx^{\prime}\\\
&-\frac{k_{\mathrm{B}}T}{D}\int^{x}_{0}\frac{\sum_{i}\beta_{i}(x^{\prime})\left(\rho_{i}(x^{\prime})Q-J_{i}\right)}{\pi
R(x^{\prime})^{2}}dx^{\prime},\end{split}$ (29)

and the pressure difference is given by $\Delta p=p(-H)-p(H)$.

##  III Perturbation theory

The applied voltage and pressure difference, $\Delta\psi$ and $\Delta p$,
generate the ionic current $I$ and the solution flux $Q$ across the
nanochannel. We consider the perturbation expansion of the fluxes as functions
of the fields given by:

$\displaystyle Q$ $\displaystyle=$ $\displaystyle D\frac{\pi
R_{0}^{2}}{L}\sum_{n+m=0}^{\infty}Q^{n}_{m}\frac{(\Delta
P)^{n}}{n!}\frac{(\Delta\varphi)^{m}}{m!},$ (30) $\displaystyle I$
$\displaystyle=$ $\displaystyle e\rho_{0}D\frac{\pi
R_{0}^{2}}{L}\sum_{n+m=0}^{\infty}I^{n}_{m}\frac{(\Delta
P)^{n}}{n!}\frac{(\Delta\varphi)^{m}}{m!},$ (31)

where $\Delta P=\Delta p/\rho_{0}k_{\mathrm{B}}T$ and
$\Delta\varphi=e\Delta\psi/k_{\mathrm{B}}T$ are the dimensionless pressure
difference and applied voltage, respectively. The summation
$\sum_{n+m=0}^{\infty}$ denotes the sum of the terms
$(n,m)=(0,0),(1,0),(0,1),(2,0),(1,1),(0,2),\cdots$. $Q_{0}^{0}$ and
$I^{0}_{0}$ are zero by definition. $Q^{1}_{0}$, $Q^{0}_{1}$, $I^{1}_{0}$, and
$I^{0}_{1}$ are the linear electrokinetic coefficients, whereas $Q^{n}_{m}$
and $I^{n}_{m}$ for $n+m\geq 2$ are the nonlinear electrokinetic coefficients.
In this study, we obtain the analytic expressions of the electrokinetic
coefficients up to the order of $n+m=3$. Detailed derivation of the framework
of the perturbation theory is explained in the following sections.

###  III.1 Derivation of expanded Poisson equations

In this section, we derive expanded Poisson equations of each order $(n,m)$.
To solve eq. 27 by the perturbation theory, it is necessary to expand eq. 27
with respect to $\mathrm{Pe}$ and $\Delta\varphi$ not $\Delta P$ and
$\Delta\varphi$. At first, we expand $\varphi(x)$, $W_{i}(x)$, and $c_{i}(x)$
with respect to $\mathrm{Pe}$ and $\Delta\varphi$,

$\displaystyle\varphi(x)$ $\displaystyle=$
$\displaystyle\sum_{n+m=0}^{\infty}{\varphi}^{n}_{m}(x)\frac{(\mathrm{Pe})^{n}}{n!}\frac{(\Delta\varphi)^{m}}{m!},$
(32) $\displaystyle W_{i}(x)$ $\displaystyle=$
$\displaystyle\sum_{n+m=0}^{\infty}{W_{i}}^{n}_{m}(x)\frac{(\mathrm{Pe})^{n}}{n!}\frac{(\Delta\varphi)^{m}}{m!},$
(33) $\displaystyle c_{i}(x)$ $\displaystyle=$
$\displaystyle\sum_{n+m=0}^{\infty}{c_{i}}^{n}_{m}(x)\frac{(\mathrm{Pe})^{n}}{n!}\frac{(\Delta\varphi)^{m}}{m!}.$
(34)

Using eq.17, the expanded ${W_{i}}^{n}_{m}(x)$ is given by:

${W_{i}}^{n}_{m}(x)=\left\\{\begin{array}[]{ll}q_{i}\varphi^{n}_{m}(x)&\textrm{for
}(n,m)\neq(1,0),\\\ \\\ \displaystyle q_{i}\varphi^{1}_{0}(x)-\frac{\pi
R_{0}^{2}}{L}&\displaystyle\int^{x}_{0}\frac{\beta_{i}(x^{\prime})}{\pi
R(x^{\prime})^{2}}dx^{\prime}\\\ \\\ &\textrm{for
}(n,m)=(1,0).\end{array}\right.$ (35)

To calculate ${c_{i}}^{n}_{m}(x)$, we need to expand ${f_{i}}(x)$ as:

$f_{i}(x)=\frac{L}{\pi{R_{0}}^{2}}\sum_{n+m=0}^{\infty}{f_{i}}^{n}_{m}(x)\frac{(\mathrm{Pe})^{n}}{n!}\frac{(\Delta\varphi)^{m}}{m!}.$
(36)

Then, we obtain the expansion for ${f_{i}}(x)/{f_{i}}(H)$,

$\frac{f_{i}(x)}{f_{i}(H)}=\frac{{f_{i}}^{0}_{0}(x)}{{f_{i}}^{0}_{0}(H)}+\frac{{f_{i}}^{1}_{0}(x)}{{f_{i}}^{0}_{0}(H)}\mathrm{Pe}+\frac{{f_{i}}^{0}_{1}(x)}{{f_{i}}^{0}_{0}(H)}\Delta\varphi+\cdots,$
(37)

where we use ${f_{i}}^{1}_{0}(H)={f_{i}}^{0}_{1}(H)=0$. Then, we obtain
${c_{i}}^{n}_{m}(x)$ up to the second order $n+m=2$,

$\displaystyle{c_{i}}^{0}_{0}(x)$ $\displaystyle=$ $\displaystyle 1,$ (38)
$\displaystyle{c_{i}}^{1}_{0}(x)$ $\displaystyle=$ $\displaystyle
p_{i}\left(\frac{1}{2}-h{f_{i}}^{0}_{0}(x)\right),$ (39)
$\displaystyle{c_{i}}^{0}_{1}(x)$ $\displaystyle=$ $\displaystyle
q_{i}\left(\frac{1}{2}-h{f_{i}}^{0}_{0}(x)\right),$ (40)
$\displaystyle{c_{i}}^{2}_{0}(x)$ $\displaystyle=$
$\displaystyle\frac{p_{i}^{2}}{4}-2p_{i}h{f_{i}}^{1}_{0}(x),$ (41)
$\displaystyle{c_{i}}^{1}_{1}(x)$ $\displaystyle=$
$\displaystyle\frac{p_{i}q_{i}}{4}-q_{i}h{f_{i}}^{1}_{0}(x)-p_{i}h{f_{i}}^{0}_{1}(x),$
(42) $\displaystyle{c_{i}}^{0}_{2}(x)$ $\displaystyle=$
$\displaystyle\frac{1}{4}-2q_{i}h{f_{i}}^{0}_{1}(x).$ (43)

where $h=1/{f_{i}}^{0}_{0}(H)$ is $i$-independent. For convenience of the
following analysis, we define a functional $F[x;g(x)]$ by:

$F[x;g(x)]=\frac{\pi{R_{0}}^{2}}{L}\int^{x}_{-H}\frac{g(x^{\prime})dx^{\prime}}{\pi
R(x^{\prime})^{2}}.$ (44)

Then, we obtain the expressions for ${f_{i}}^{n}_{m}(x)$ up to the second
order given by:

$\displaystyle{f_{i}}^{0}_{0}(x)$ $\displaystyle=$ $\displaystyle F(x;1),$
(45) $\displaystyle{f_{i}}^{1}_{0}(x)$ $\displaystyle=$ $\displaystyle
F(x;{W_{i}}^{1}_{0}(x)),$ (46) $\displaystyle{f_{i}}^{0}_{1}(x)$
$\displaystyle=$ $\displaystyle F(x;{W_{i}}^{0}_{1}(x)),$ (47)
$\displaystyle{f_{i}}^{2}_{0}(x)$ $\displaystyle=$ $\displaystyle
F(x;{W_{i}}^{2}_{0}(x)+{W_{i}}^{1}_{0}(x)^{2}),$ (48)
$\displaystyle{f_{i}}^{1}_{1}(x)$ $\displaystyle=$ $\displaystyle
F(x;{W_{i}}^{1}_{1}(x)+{W_{i}}^{1}_{0}(x){W_{i}}^{0}_{1}(x)),$ (49)
$\displaystyle{f_{i}}^{0}_{2}(x)$ $\displaystyle=$ $\displaystyle
F(x;{W_{i}}^{0}_{2}(x)+{W_{i}}^{0}_{1}(x)^{2}),$ (50)

where we use ${W_{i}}^{0}_{0}(x)=0$. For the zeroth order,
${c_{i}}^{0}_{0}(x)=1$ and ${W_{i}}^{0}_{0}(x)=q_{i}\varphi^{0}_{0}(x)$ give
the zeroth-order Poisson equation (usual Poisson-Boltzmann equation),

$\frac{1}{\kappa^{2}R^{2}(x)}\frac{d}{dx}\left[R^{2}(x)\frac{d\varphi^{0}_{0}(x)}{dx}\right]=\sinh\varphi^{0}_{0},$
(51)

with the boundary condition $\varphi^{0}_{0}(-H)=\varphi^{0}_{0}(H)=0$. The
solution of eq. 51 is $\varphi^{0}_{0}(x)=0$. The higher order expanded
Poisson equations are given by:

$\frac{1}{\kappa^{2}R^{2}(x)}\frac{d}{dx}\left[R^{2}(x)\frac{d\varphi^{n}_{m}(x)}{dx}\right]=-\rho^{n}_{m}(x),$
(52)

where $\rho^{n}_{m}(x)=\frac{1}{2}\sum_{i}q_{i}{\rho_{i}}^{n}_{m}(x)$ is
induced charge density, and the expressions for ${\rho_{i}}^{n}_{m}(x)$ up to
the second order are given by:

$\displaystyle{\rho_{i}}^{1}_{0}(x)$ $\displaystyle=$
$\displaystyle{c_{i}}^{1}_{0}(x)-{W_{i}}^{1}_{0}(x),$ (53)
$\displaystyle{\rho_{i}}^{0}_{1}(x)$ $\displaystyle=$
$\displaystyle{c_{i}}^{0}_{1}(x)-{W_{i}}^{0}_{1}(x),$ (54)
$\displaystyle{\rho_{i}}^{2}_{0}(x)$ $\displaystyle=$
$\displaystyle{c_{i}}^{2}_{0}(x)-2{c_{i}}^{1}_{0}(x){W_{i}}^{1}_{0}(x)-{W_{i}}^{2}_{0}(x)+{W_{i}}^{1}_{0}(x)^{2},$
$\displaystyle{\rho_{i}}^{1}_{1}(x)$ $\displaystyle=$
$\displaystyle{c_{i}}^{1}_{1}(x)-{c_{i}}^{1}_{0}(x){W_{i}}^{0}_{1}(x)-{c_{i}}^{0}_{1}(x){W_{i}}^{1}_{0}(x)$
(56) $\displaystyle-{W_{i}}^{1}_{1}(x)+{W_{i}}^{0}_{1}(x){W_{i}}^{1}_{0}(x),$
$\displaystyle{\rho_{i}}^{0}_{2}(x)$ $\displaystyle=$
$\displaystyle{c_{i}}^{0}_{2}(x)-2{c_{i}}^{0}_{1}(x){W_{i}}^{0}_{1}(x)-{W_{i}}^{0}_{2}(x)+{W_{i}}^{0}_{1}(x)^{2}.$

In the derivation of eqs. 53 to III.1, we use ${c_{i}}^{0}_{0}(x)=1$,
${W_{i}}^{0}_{0}(x)=0$, and

$\begin{split}\mathrm{e}^{\pm
W_{i}(x)}&=1\pm{W_{i}}^{1}_{0}(x)\mathrm{Pe}\pm{W_{i}}^{0}_{1}(x)\Delta\varphi\\\
&+[{W_{i}}^{1}_{0}(x)^{2}\pm{W_{i}}^{2}_{0}(x)]\frac{(\mathrm{Pe})^{2}}{2}\\\
&+[{W_{i}}^{1}_{0}(x){W_{i}}^{0}_{1}(x)\pm{W_{i}}^{1}_{1}(x)](\mathrm{Pe})(\Delta\varphi)\\\
&+[{W_{i}}^{0}_{1}(x)^{2}\pm{W_{i}}^{0}_{2}(x)]\frac{(\Delta\varphi)^{2}}{2}+\cdots.\end{split}$
(58)

The boundary conditions of eq. 52 are
$\varphi^{n}_{m}(-H)=\varphi^{n}_{m}(H)=0$ for $(n,m)\neq(0,1)$ and
$\varphi^{0}_{1}(-H)=-\varphi^{0}_{1}(H)=1/2$ for $(n,m)=(0,1)$. Because
$\rho^{n}_{m}(x)$ ($n+m\geq 1$) has a linear term $-\varphi^{n}_{m}(x)$, we
can extract the source term of the linear differential equations by
$S^{n}_{m}(x)=\rho^{n}_{m}(x)+\varphi^{n}_{m}(x)$ which are given by:

$\displaystyle S^{1}_{0}(x)=\frac{1}{2}$
$\displaystyle\times\left[{c_{+}}^{1}_{0}(x)-{c_{-}}^{1}_{0}(x)+\frac{\pi{R_{0}}^{2}}{L}\int^{x}_{0}\frac{\beta_{+}(x^{\prime})-\beta_{-}(x^{\prime})}{\pi
R(x^{\prime})^{2}}dx^{\prime}\right],$ (59) $\displaystyle
S^{0}_{1}(x)=\frac{1}{2}\left[{c_{+}}^{0}_{1}(x)-{c_{-}}^{0}_{1}(x)\right],$
(60) $\displaystyle
S^{2}_{0}(x)=\sum_{i=\pm}\frac{q_{i}}{2}\left[{c_{i}}^{2}_{0}(x)+{W_{i}}^{1}_{0}(x)^{2}-2{c_{i}}^{1}_{0}(x){W_{i}}^{1}_{0}(x)\right],$
(61) $\displaystyle
S^{1}_{1}(x)=\sum_{i=\pm}\frac{q_{i}}{2}\left[{c_{i}}^{1}_{1}(x)+{W_{i}}^{1}_{0}(x){W_{i}}^{0}_{1}(x)\right.$
$\displaystyle\qquad\qquad\qquad\qquad\left.-{c_{i}}^{1}_{0}(x){W_{i}}^{0}_{1}(x)-{c_{i}}^{0}_{1}(x){W_{i}}^{1}_{0}(x)\right],$
(62) $\displaystyle
S^{0}_{2}(x)=\sum_{i=\pm}\frac{q_{i}}{2}\left[{c_{i}}^{0}_{2}(x)+{W_{i}}^{0}_{1}(x)^{2}-2{c_{i}}^{0}_{1}(x){W_{i}}^{0}_{1}(x)\right].$

###  III.2 Solution of expanded Poisson equation in each order

The first-order and the second-order Poisson equations are obtained as:

$\frac{1}{\kappa^{2}R^{2}(x)}\frac{d}{dx}\left[R^{2}(x)\frac{d\varphi^{n}_{m}(x)}{dx}\right]-\varphi^{n}_{m}(x)=-S^{n}_{m}(x),$
(64)

Further calculation yields analytical expressions for $S^{1}_{0}(x)$ and
$S^{0}_{1}(x)$ as follows:

$\displaystyle S^{1}_{0}(x)$ $\displaystyle=$
$\displaystyle\left\\{\begin{array}[]{ll}\displaystyle-\mathrm{sgn}(x)\left[-\frac{\sigma_{0}}{2}h\nu
t(x)\right],&\displaystyle\textrm{for }|x|>\frac{L}{2}\\\ \\\
\displaystyle\sigma_{0}h\nu t^{\prime}\frac{x}{L},&\displaystyle\textrm{for
}|x|<\frac{L}{2}\end{array},\right.$ (68) $\displaystyle S^{0}_{1}(x)$
$\displaystyle=$
$\displaystyle\left\\{\begin{array}[]{ll}\displaystyle-\mathrm{sgn}(x)\left[\frac{1}{2}-h\nu
t(x)\right],&\displaystyle\textrm{for }|x|>\frac{L}{2}\\\ \\\
\displaystyle-h\frac{x}{L},&\displaystyle\textrm{for
}|x|<\frac{L}{2}\end{array},\right.$ (72)

where $\sigma_{0}=\sum_{i}q_{i}\alpha_{i}$, $\nu=\delta/L$,
$t(x)=\mathrm{e}^{-\frac{|x|-L/2}{\delta}}-\mathrm{e}^{\frac{H-L/2}{\delta}}$,
and $t^{\prime}=t(\frac{L}{2})$. Using this notation, we have $h=1/(1+2\nu
t^{\prime})$. For $(n,m)=(1,0)$ and $(0,1)$ we obtain the approximated
solutions,

$\displaystyle\varphi^{1}_{0}(x)$ $\displaystyle=$ $\displaystyle
S^{1}_{0}(x),$ (73) $\displaystyle\varphi^{0}_{1}(x)$ $\displaystyle=$
$\displaystyle S^{0}_{1}(x).$ (74)

The exact solution is given by supplementary material.

Using the approximate solutions for $\varphi^{1}_{0}(x)$ and
$\varphi^{0}_{1}(x)$, the analytical expressions for $S^{2}_{0}(x)$,
$S^{1}_{1}(x)$, and $S^{0}_{2}(x)$ are obtained as follows:

$\displaystyle S^{2}_{0}(x)$ $\displaystyle=$
$\displaystyle\left\\{\begin{array}[]{ll}\displaystyle\frac{\sigma_{0}(2-\sigma_{1})}{4}h^{2}\nu^{2}t(x)^{2},&\displaystyle\textrm{for
}|x|>\frac{L}{2},\\\ \\\ \displaystyle\frac{\sigma_{0}(2-\sigma_{1})}{8}h\nu
t^{\prime}-\frac{\sigma_{0}(2-\sigma_{1})}{2}&\displaystyle h^{2}\nu
t^{\prime}\frac{x^{2}}{L^{2}},\\\ &\displaystyle\textrm{for
}|x|<\frac{L}{2},\end{array}\right.$ (79) $\displaystyle S^{1}_{1}(x)$
$\displaystyle=$
$\displaystyle\left\\{\begin{array}[]{ll}\displaystyle\frac{2-\sigma_{1}}{4}h^{2}\nu^{2}t(x)^{2},&\displaystyle\textrm{for
}|x|>\frac{L}{2},\\\ \\\ \displaystyle\frac{2-\sigma_{1}}{8}h\nu
t^{\prime}-\frac{2-\sigma_{1}}{2}h^{2}\nu
t^{\prime}\frac{x^{2}}{L^{2}},&\displaystyle\textrm{for
}|x|<\frac{L}{2},\end{array}\right.$ (83) $\displaystyle S^{0}_{2}(x)$
$\displaystyle=$ $\displaystyle 0,$ (84)

where $\sigma_{1}=\sum_{i}\alpha_{i}$. For $(n,m)=(2,0)$,$(1,1)$, and $(0,2)$,
we obtain the approximated solutions,

$\displaystyle\varphi^{2}_{0}(x)$ $\displaystyle=$ $\displaystyle
S^{2}_{0}(x),$ (85) $\displaystyle\varphi^{1}_{1}(x)$ $\displaystyle=$
$\displaystyle S^{1}_{1}(x),$ (86) $\displaystyle\varphi^{0}_{2}(x)$
$\displaystyle=$ $\displaystyle 0.$ (87)

Derivation of the exact solution is described in supplementary material. Using
eqs. 68–87, the induced charge density $\rho^{n}_{m}(x)$ of each order inside
the nanochannel is obtained as:

$\displaystyle\rho^{1}_{0}(x)$ $\displaystyle=$ $\displaystyle 0,$ (88)
$\displaystyle\rho^{0}_{1}(x)$ $\displaystyle=$ $\displaystyle 0,$ (89)
$\displaystyle\rho^{2}_{0}(x)$ $\displaystyle=$
$\displaystyle\frac{\sigma_{0}(2-\sigma_{1})}{\kappa^{2}}h^{2}\nu
t^{\prime}\frac{1}{L^{2}},$ (90) $\displaystyle\rho^{1}_{1}(x)$
$\displaystyle=$ $\displaystyle\frac{2-\sigma_{1}}{\kappa^{2}}h^{2}\nu
t^{\prime}\frac{1}{L^{2}},$ (91) $\displaystyle\rho^{0}_{2}(x)$
$\displaystyle=$ $\displaystyle 0.$ (92)

After obtaining $\varphi^{n}_{m}(x)$ and $\rho^{n}_{m}(x)$ up to the second
order, we can calculate $\Delta p^{n}_{m}$ and ${J_{i}}^{n}_{m}$ up to the
third order. The individual ionic current $J_{i}$ defined in eq. 26 is
expanded with respect to $\mathrm{Pe}$ and $\Delta\varphi$ as:

$J_{i}(\mathrm{Pe},\Delta\varphi)=\rho_{0}D\frac{\pi
R_{0}^{2}}{L}\sum_{n+m=0}^{\infty}{J_{i}}^{n}_{m}\frac{(\mathrm{Pe})^{n}}{n!}\frac{(\Delta\varphi)^{m}}{m!}.$
(93)

We also calculate $\Delta p^{n}_{m}$ which is defined by:

$\Delta
p(\mathrm{Pe},\Delta\varphi)=\rho_{0}k_{\mathrm{B}}T\sum_{n+m=0}^{\infty}\Delta
p^{n}_{m}\frac{(\mathrm{Pe})^{n}}{n!}\frac{(\Delta\varphi)^{m}}{m!}.$ (94)

Explicit formulae for ${J_{i}}^{n}_{m}$ and $\Delta p^{n}_{m}$ are available
in supplementary material.

Because our framework of the perturbation theory is based on the expansion
with $\mathrm{Pe}$ and $\Delta\varphi$, we need to change the variables into
$\Delta P$ and $\Delta\varphi$. The detailed procedure and explicit formulae
for $Q^{n}_{m}$ and $I^{n}_{m}$ are given in supplementary material.

##  IV Results and Discussion

In this section, we investigate typical parameter sets of nanochannels and
compare full numerical solutions with analytical results of the perturbation
theory. Fig. 2 shows the voltage-dependent streaming current
$I_{\mathrm{str}}(\Delta p,\Delta\psi)=I(\Delta p,\Delta\psi)-I(0,\Delta\psi)$
which was recently measured in the experiments. [7, 10] The ion-wall frictions
are set at $\alpha_{+}=0.002$ and $\alpha_{-}=0.001$. Other parameters are set
at $L=1\,\mu$m, $R_{0}=1\,$nm, $D=2\times 10^{-9}\,$m2/s, $\eta=1\,$mPa/s,
$T=300\,$K, $\varepsilon=80$, $\rho_{0}=1\,$M, $\delta=250\,$nm, and
$2H=1.25\,\mu$m. The lines are calculated analytically up to the terms of
$n+m=3$ of eq. 31, and the points are obtained by numerically solving the
electro-hydrodynamic equations. The different color represents the different
voltage $\Delta\psi=0$, $25$, $50$, $75$, and $100\,$mV. The analytical lines
perfectly agree with the points obtained numerically, indicating that our
analytic theory is correct. In the low hydrodynamic permeability case (a)
$b=1\,$nm, the nonlinear streaming current is not so drastic, whereas in the
high hydrodynamic permeability case (b) $b=10\,\mu$m, the streaming current
exhibits a quadratic dependence on the applied pressure difference. Because
the numerical solution includes higher-order terms $n+m>3$, these terms are
found to be negligible in this condition.

Figure 2:  The voltage-dependent streaming currents $I_{\mathrm{str}}(\Delta
p,\Delta\psi)$ are plotted for (a) low hydrodynamic permeability $b=1\,$nm and
(b) high hydrodynamic permeability $b=10\,\mu$m. The points are numerically
obtained by solving eq. 27, whereas the lines are eq. 31 with analytically
obtained coefficients up to $n+m=3$. The different color represents
$\Delta\psi=0$, $25$, $50$, $75$, and $100\,$mV.

Although the analytical expressions for all of the electrokinetic coefficients
up to the third order are complicated (eqs. S39–S50 in supplementary
material), we take a long-channel limit, $\nu=\delta/L\to 0$, and the infinite
system size ($2H\to\infty$) to derive simplified analytic expressions. Leading
terms of the linear coefficients, $Q^{1}_{0}$, $Q^{0}_{1}$, $I^{1}_{0}$, and
$I^{0}_{1}$, are obtained as follows:

$\displaystyle Q^{1}_{0}$ $\displaystyle=$ $\displaystyle K,$ (95)
$\displaystyle Q^{0}_{1}$ $\displaystyle=$ $\displaystyle
I^{1}_{0}=K(\alpha_{+}-\alpha_{-}),$ (96) $\displaystyle I^{0}_{1}$
$\displaystyle=$ $\displaystyle 2+K(\alpha_{+}-\alpha_{-})^{2},$ (97)

where $K(>0)$ is the non-dimensional hydrodynamic permeation coefficient,
which is given by:

$\frac{1}{K}=\mathrm{Pe}_{\mathrm{osm}}^{-1}+\sum_{i}\alpha_{i}(1-\alpha_{i}),$
(98)

and $\mathrm{Pe}_{\mathrm{osm}}=\rho_{0}k_{\mathrm{B}}TR_{0}b/2\eta D$ denotes
the hydrodynamic permeation coefficient without the effect of ion-wall
friction. The ion-wall friction modifies the hydrodynamic permeation
coefficient from $\mathrm{Pe}_{\mathrm{osm}}$ to $K$. Similarly, the ionic
conductance $I^{0}_{1}$ is modified from $2$ to
$2+K(\alpha_{+}-\alpha_{-})^{2}$ due to the ion-wall friction. Because
$Q^{1}_{0}$ and $I^{0}_{1}$ are diagonal components of the macroscopic
transport matrix, they are never negative as a result of positive entropy
production. The coefficients of the electro-osmosis and streaming current are
the same, $Q^{0}_{1}=I^{1}_{0}$, reflecting the Onsager reciprocal relation.
Interestingly, even though the nanochannel is uncharged in this model, the
electro-osmotic mobility (or the zeta potential) is non-zero. When the anions
are more bounded to the nanochannel than the cations, i.e.
$\alpha_{+}>\alpha_{-}$, it is intuitive that the zeta potential becomes
negative due to the difference between the ion-wall frictions. In experiments,
the nanochannel/water interface is considered to be spontaneously charged due
to the friction asymmetry as well as the ion adsorption.

All the second-order coefficients, $Q^{2}_{0}$, $Q^{1}_{1}$, $Q^{0}_{2}$,
$I^{2}_{0}$, $I^{1}_{1}$, and $I^{0}_{2}$, are strictly zero because of the
symmetric geometry of the nanochannel. All the third-order coefficients are
proportional to $\delta/L$ in the lowest order. Neglecting the higher order
terms $(\delta/L)^{2}$, the analytical expressions for the water flux are
obtained as follows:

$\displaystyle Q^{3}_{0}$ $\displaystyle=$
$\displaystyle-\frac{K^{4}}{4}\sigma_{1}\sigma_{2}(2-\sigma_{1})\frac{\delta}{L},$
(99) $\displaystyle Q^{2}_{1}$ $\displaystyle=$
$\displaystyle\frac{K^{3}}{12}\sigma_{0}\sigma_{1}(2-\sigma_{1})\left(1-3K\sigma_{2}\right)\frac{\delta}{L},$
(100) $\displaystyle Q^{1}_{2}$ $\displaystyle=$
$\displaystyle\frac{K^{3}}{12}{\sigma_{0}}^{2}\sigma_{1}(2-\sigma_{1})\left(2-3K\sigma_{2}\right)\frac{\delta}{L},$
(101) $\displaystyle Q^{0}_{3}$ $\displaystyle=$
$\displaystyle\frac{K^{3}}{4}{\sigma_{0}}^{3}\sigma_{1}(2-\sigma_{1})\left(1-K\sigma_{2}\right)\frac{\delta}{L},$
(102)

where $\sigma_{2}=\sum_{i}\alpha_{i}(1-\alpha_{i})$, and for the ionic
current,

$\displaystyle I^{3}_{0}$ $\displaystyle=$
$\displaystyle\frac{K^{3}}{4}\sigma_{0}\sigma_{1}(2-\sigma_{1})\left(1-K\sigma_{2}\right)\frac{\delta}{L},$
(103) $\displaystyle I^{2}_{1}$ $\displaystyle=$
$\displaystyle\frac{K^{2}}{12}\sigma_{1}(2-\sigma_{1})\left[2+K{\sigma_{0}}^{2}(4-3K\sigma_{2})\right]\frac{\delta}{L},$
(104) $\displaystyle I^{1}_{2}$ $\displaystyle=$
$\displaystyle\frac{K^{2}}{12}\sigma_{0}\sigma_{1}(2-\sigma_{1})\left[4+K{\sigma_{0}}^{2}(5-3K\sigma_{2})\right]\frac{\delta}{L},$
(105) $\displaystyle I^{0}_{3}$ $\displaystyle=$
$\displaystyle\frac{K^{2}}{4}{\sigma_{0}}^{2}\sigma_{1}(2-\sigma_{1})\left[2+K{\sigma_{0}}^{2}(2-K\sigma_{2})\right]\frac{\delta}{L}.$
(106)

Eqs. 99 to 106 show that all the third order coefficients are proportional to
$2-\sigma_{1}$, and $K^{2}$. This clearly demonstrates that the large
hydrodynamic permeation, $K\gg 1$, is necessary for strong third-order
responses. Furthermore, zero ion-wall friction
($\sigma_{1}=\alpha_{+}+\alpha_{-}=2$) remove all the third-order responses.
According to the Onsager reciprocal relation for nonlinear transport, [22] the
relations: $Q^{2}_{1}=I^{3}_{0}$, $Q^{1}_{2}=I^{2}_{1}$, and
$Q^{0}_{3}=I^{1}_{2}$, are expected. However, these relations do not hold in
eqs. 99 to 106 because our framework includes the responses of the induced
charge density for $n+m=2$, which is discussed in the next paragraph. In fact,
these relation did not hold in the previous membrane experiments, [22] and
thus, we think our theory correctly reflects the properties of nonlinear
responses.

Figure 3:  Analytically calculated $\varphi^{n}_{m}(x)$ and $\rho^{n}_{m}(x)$
for $n+m=1$ and $2$. (a) The profile of $\varphi^{1}_{0}(x)$, (b)
$\varphi^{0}_{1}(x)$, (c) $\varphi^{0}_{2}(x)$, (d) $\varphi^{2}_{0}(x)$ and
$\varphi^{1}_{1}(x)$, (e) $\rho^{1}_{0}(x)$, $\rho^{0}_{1}(x)$, and
$\rho^{0}_{2}(x)$, (f) $\rho^{2}_{0}(x)$ and $\rho^{1}_{1}(x)$. The vertical
axes are arbitrary unit except for $\varphi^{0}_{1}(x)$. For $\rho^{2}_{0}(x)$
and $\rho^{1}_{1}(x)$, we plot
$R(x)^{-2}\frac{d}{dx}\left[R(x)^{2}\frac{d}{dx}\right]S^{2}_{0}(x)$ because
they are proportional to it. The used parameters are the same as Fig. 2.

In Fig. 3, we plot analytically calculated $\varphi^{n}_{m}(x)$ and
$\rho^{n}_{m}(x)$. The first order term $\varphi^{1}_{0}(x)$ and
$\varphi^{0}_{1}(x)$ exhibit linear slopes inside the nanochannel (Fig. 3ab),
indicating that the linear response does not affect the induced charge,
$\rho^{n}_{m}(x)$, (Fig. 3e). The $(\Delta\varphi)^{2}$ term of $\varphi(x)$,
$\varphi^{0}_{2}(x)$, is zero (Fig. 3c), indicating that the external electric
field without the flux does not induce the charge inside the nanochannel.
However, the other second-order potentials, $\varphi^{2}_{0}(x)$ and
$\varphi^{1}_{1}(x)$, exhibit parabolas inside the nanochannel (Fig. 3d),
demonstrating that the external flow $Q$ induces the uniform electrification
inside the nanochannel (Fig. 3f) because $\mathrm{Pe}$ is linear with $Q$ (see
eq. S33). Although the sign and the magnitude of the induced charge due to the
external flow $Q$ are complicated, the induced charge inside the nanochannel
is proportional to $2-\sigma_{1}$ (see eqs. 90 and 91). This fact provides
intuitive understanding of the nonlinear electrokinetics: the third-order
electrokinetics in the nanochannel is attributed to the uniform induced charge
proportional to $\mathrm{Pe}^{2}$ and $\mathrm{Pe}\cdot\Delta\varphi$. In
terms of the external fields $\Delta P$ and $\Delta\varphi$, the induced
charge inside the nanochannel is proportional to $(\Delta
P)^{n}(\Delta\varphi)^{m}$ with $n+m=2$ because $Q$ is linearly coupled with
$\Delta\varphi$ and $\Delta P$ in the lowest order. As a result, the
electrokinetic flows, $I$ and $Q$, have a nonlinear term $(\Delta
P)^{n}(\Delta\varphi)^{m}$ with $n+m=3$.

Table 1:  Linear and nonlinear electrokinetic coefficients extracted from the previous nanochannel experiments [10]. The electrokinetic coefficients, $I^{n}_{m}$, are converted from the experimental data. The numbers refer those described in Supplementary Information of Ref. 10. The salt concentration is $\rho_{0}=1\,$M and $\mathrm{pH}=5.5$. Extraction detail is described in supplementary material.  | $I^{1}_{0}$ | $I^{0}_{1}$ | $I^{3}_{0}$ | $I^{2}_{1}$ | $I^{1}_{2}$ | $I^{0}_{3}$
---|---|---|---|---|---|---
#2 |  $\sim 0$111negligibly small compared to $I^{2}_{1}$ | 48 | $\sim 0^{\mathrm{a}}$ | 4000 | $\sim 0^{\mathrm{a}}$ |  $\sim 0$222negligibly small compared to $I^{0}_{1}$
#3 |  $\sim 0$a | 32 | $\sim 0^{\mathrm{a}}$ | 93000 | $\sim 0^{\mathrm{a}}$ | $\sim 0^{\mathrm{b}}$

Discussion about relevant experiments [10] is crucially important for the
theoretical work. In Ref. 10, the measurements of the streaming current were
performed in a single carbon nanotube with a radius $2\,$nm and a length about
$1\,\mu$m, [10] where the ionic current exhibits a quadratic dependence on the
pressure difference, $(\Delta p)^{2}\Delta\psi$. In their paper, this
quadratic dependence was qualitatively explained by the similar electro-
hydrodynamic model with the energy-barrier effect due to the confinement. [10]
Here, we quantitatively discuss this experiment by the ion-wall friction.

Table. 1 lists the electrokinetic coefficients, $I^{n}_{m}$, extracted from
Ref. 10. The ionic current of the single carbon nanotube with varied
$\Delta\psi$ and $\Delta p=0$ was perfectly linear with $\Delta\psi$ in the
range of $-50\,$mV to $50\,$mV, suggesting $I^{0}_{3}$ is negligible compared
to $I^{0}_{1}$. Furthermore, the voltage-dependent streaming current,
$I_{\mathrm{str}}(\Delta p,\Delta\psi)=I(\Delta p,\Delta\psi)-I(0,\Delta\psi)$
exhibits the parabolic dependence on $\Delta p$, and the linear streaming
current was negligible for finite $\Delta\psi$. Therefore, $I^{1}_{0}$,
$I^{1}_{2}$, $I^{3}_{0}$ are negligible compared to $I^{2}_{1}$. To extract
the friction parameters $K$ and $\alpha_{i}$ from $I^{1}_{2}$, it is clear
that the characteristic length at the joint between the nanochannel and
reservoirs is significantly responsible for the determination of the nonlinear
electrokinetics. The importance of the joint length $\delta$ has not been
recognized so far, and how long the effective $\delta$ is in the actual
experiment remains unknown. The obtained dimensionless conductance,
$I^{0}_{1}$, is in the order of $10$, suggesting that our theory yields $K\sim
10^{2}$ from eq. 97 because $|\sigma_{0}|=|\alpha_{+}-\alpha_{-}|<1$.
Furthermore, the nonlinear coefficient $I^{2}_{1}$ is in the order of $10^{3}$
to $10^{5}$, which is a significantly strong nonlinear response. Using eq.
104, this is achieved only when $K\gg 1$ and $\sigma_{1}\neq 0$ nor $2$. A
complete fit of all the friction parameters in the model with the experimental
data is exciting, and more experiments on the water flow $Q$ are necessary.

##  V Summary

We have constructed the perturbation theory to calculate the nonlinear
electrokinetic coefficients of confined electrolyte solution in a nanochannel.
The theory is applied to the model electro-hydrodynamic equations, which
include the water-wall and ion-wall friction effects. The linear and nonlinear
electrokinetic coefficients up to the third order are analytically calculated,
and simple analytical expressions for the electrokinetic coefficients are
obtained in the long channel limit. Furthermore, the responses of the
electrostatic potential profiles along the nanochannel are analytically
obtained, which reveals that the uniform flow-induced charge inside the
nanochannel is attributed to the ion-wall friction. As a result, the induced
charge inside the channel has a quadratic dependence on the external fields,
and thus, the third-order nonlinear electrokinetics emerges.

Our framework is general, and thus, it is easy to include non-uniform
diffusion profiles under confinement. In this study, additional relevant
effects, namely, the asymmetric geometry of the nanochannel, the fixed surface
charge on the nanochannel surface, and the energy barrier of ions inside the
nanochannel, are not considered within the model. Further inclusion of these
effects and quantitative comparison of the water flows and ionic currents
between theory and experiments could be a fruitful frontier of future
nanofluidics research.

## SUPPLEMENTARY MATERIAL

See supplementary material for derivation of electrohydrodynamic equations via
Onsager variational principle, the exact solution of expanded Poisson
equations, calculation of $\Delta p^{n}_{m}$ and ${J_{i}}^{n}_{m}$ up to
$n+m=3$ and conversion to $Q^{n}_{m}$ and $I^{n}_{m}$, and analysis on the
previous experimental data.

###### Acknowledgements.

This work was supported by JST, Presto Grant number JPMJPR21O2, JSPS
Fellowships 201860001, and JSPS Grant 18KK0151, 20K14430. The author thanks
Lydéric Bocquet for the host of his stay in Paris in 2019 where he started
this work and acceptance to submit it alone.

## References

  * Kim  _et al._ [2007] S. J. Kim, Y.-C. Wang, J. H. Lee, H. Jang, and J. Han, “Concentration polarization and nonlinear electrokinetic flow near a nanofluidic channel,” Phys. Rev. Lett. 99, 044501 (2007).
  * Bocquet and Charlaix [2010] L. Bocquet and E. Charlaix, “Nanofluidics, from bulk to interfaces,” Chem. Soc. Rev. 39, 1073–1095 (2010).
  * Majumder  _et al._ [2005] M. Majumder, N. Chopra, R. Andrews, and B. J. Hinds, “Enhanced flow in carbon nanotubes,” Nature 438, 44 (2005).
  * Holt  _et al._ [2006] J. K. Holt, H. G. Park, Y. Wang, M. Stadermann, A. B. Artyukhin, C. P. Grigoropoulos, A. Noy, and O. Bakajin, “Fast mass transport through sub-2-nanometer carbon nanotubes,” Science 312, 1034–1037 (2006).
  * Secchi  _et al._ [2016a] E. Secchi, S. Marbach, A. Niguès, D. Stein, A. Siria, and L. Bocquet, “Massive radius-dependent flow slippage in carbon nanotubes,” Nature 537, 210–213 (2016a).
  * Esfandiar  _et al._ [2017] 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).
  * Mouterde  _et al._ [2019] T. Mouterde, A. Keerthi, A. Poggioli, S. Dar, A. Siria, A. Geim, L. Bocquet, and B. Radha, Nature 567, 87 (2019).
  * Secchi  _et al._ [2016b] E. Secchi, A. Nigues̀, L. Jubin, A. Siria, and L. Bocquet, “Scaling behavior for ionic transport and its fluctuations in individual carbon nanotubes,” Phys. Rev. Lett. 116, 154501 (2016b).
  * Poggioli, Siria, and Bocquet [2019] A. R. Poggioli, A. Siria, and L. Bocquet, “Beyond the tradeoff: Dynamic selectivity in ionic transport and current rectification,” J. Phys. Chem. B 123, 1171–1185 (2019).
  * Marcotte  _et al._ [2020] A. Marcotte, T. Mouterde, A. Niguès, A. Siria, and L. Bocquet, “Mechanically activated ionic transport across single-digit carbon nanotubes,” Nature Materials 19, 1057–1061 (2020).
  * Davis  _et al._ [2020] S. J. Davis, M. Macha, A. Chernev, D. M. Huang, A. Radenovic, and S. Marion, “Pressure-induced enlargement and ionic current rectification in symmetric nanopores,” Nano Lett. 20, 8089–8095 (2020).
  * Yao  _et al._ [2021] Y.-C. Yao, Z. Li, A. J. Gillen, S. Yosinski, M. A. Reed, and A. Noy, “Electrostatic gating of ion transport in carbon nanotube porins: A modeling study,” J. Chem. Phys. 154, 204704 (2021).
  * Koyama  _et al._ [2021] S. Koyama, D. Inoue, A. Okada, and H. Yoshida, “Electro-osmotic diode based on colloidal nano-valves between double membranes,” Phys. Rev. Research 3 (2021).
  * Jiang  _et al._ [2022] X. Jiang, C. Zhao, Y. Noh, Y. Xu, Y. Chen, F. Chen, L. Ma, W. Ren, N. R. Aluru, and J. Feng, “Nonlinear electrohydrodynamic ion transport in graphene nanopores,” Science Advances 8, eabj2510 (2022).
  * Boon  _et al._ [2022] W. Q. Boon, T. E. Veenstra, M. Dijkstra, and R. van Roij, “Pressure-sensitive ion conduction in a conical channel: Optimal pressure and geometry,” Physics of Fluids 34, 101701 (2022).
  * You  _et al._ [2022] Y. You, A. Ismail, G.-H. Nam, S. Goutham, A. Keerthi, and B. Radha, “Angstrofluidics: Walking to the limit,” Annu. Rev. Mater. Res. 52, 189–218 (2022).
  * Bazant and Squires [2004] M. Z. Bazant and T. M. Squires, Phys. Rev. Lett. 92, 066101 (2004).
  * Tottori  _et al._ [2019] S. Tottori, K. Misiunas, U. F. Keyser, and D. J. Bonthuis, Phys. Rev. Lett. 123, 014502 (2019).
  * Khair [2022] A. S. Khair, “Nonlinear electrophoresis of colloidal particles,” Current Opinion in Colloid and Interface Science 59, 101587 (2022).
  * Squires and Bazant [2006] T. M. Squires and M. Z. Bazant, “Breaking symmetries in induced-charge electro-osmosis and electrophoresis,” J. Fluid Mech. 560, 65–101 (2006).
  * Schnitzer and Yariv [2014] O. Schnitzer and E. Yariv, “Nonlinear electrophoresis at arbitrary field strengths: small-dukhin-number analysis,” Phys. Fluids 26, 122002 (2014).
  * Rastogi, Srivastava, and Singh [1993] R. P. Rastogi, R. C. Srivastava, and S. N. Singh, Chem. Rev. 93, 1945–1990 (1993).
  * Kilic, Bazant, and Ajdari [2007] M. S. Kilic, M. Z. Bazant, and A. Ajdari, “Steric effects in the dynamics of electrolytes at large applied voltages. I. Double-layer charging,” Phys. Rev. E 75 (2007).
  * Loche  _et al._ [2019] P. Loche, C. Ayaz, A. Schlaich, Y. Uematsu, and R. R. Netz, “Giant axial dielectric response in water-filled nanotubes and effective electrostatic ion–ion interactions from a tensorial dielectric model,” J. Phys. Chem. B 123, 10850–10857 (2019).
  * Fornasiero  _et al._ [2008] F. Fornasiero, H. G. Park, J. K. Holt, M. Stadermann, C. P. Grigoropoulos, A. Noy, and O. Bakajin, “Ion exclusion by sub-2-nm carbon nanotube pores,” Proc. Natl. Acad. Sci. U.S.A. 105, 17250–17255 (2008).
  * Doi [2011] M. Doi, “Onsager’s variational principle in soft matter,” J. Phys.: Condens. Matter 23, 284118 (2011).

Supporting Information: Analytic theory of nonlinearly coupled electrokinetics
in nanochannels

##  VI Derivation of electrohydrodynamic equations via Onsager variational
principle

Eqs. 3 and 4 in the main article can be derived by the Onsager variational
principle [1] by minimizing a functional
$\mathcal{R}[v_{i}(x),v_{\mathrm{w}}(x)]$ with respect to $v_{i}(x)$ and
$v_{\mathrm{w}}(x)$. The function $\mathcal{R}$ is defined by

$\begin{split}\mathcal{R}[v_{i}(x),v_{\mathrm{w}}(x)]&=\frac{1}{2}\int\left[\sum_{i}\rho_{i}\left[\xi_{i}(v_{i}-v_{\mathrm{w}})^{2}+\lambda_{i}v_{i}^{2}\right]\pi
R(x)^{2}+\lambda_{0}v_{\mathrm{w}}^{2}2\pi R(x)\right]dx\\\
&+\frac{d}{dt}\int\sum_{i}\rho_{i}\mu_{i}\pi R(x)^{2}dx-\int
p(x)\frac{d}{dx}\left[\pi R(x)^{2}v_{\mathrm{w}}\right]dx,\end{split}$ (S1)

where the first term is the dissipation function, the second term is the time
derivative of the free energy, and the last term is the constraint of the
incompressible condition, eq. 1.

##  VII The exact solution of expanded Poisson equations

The expanded Poisson equations can be written in the form of

$\frac{1}{\kappa^{2}R(x)^{2}}\frac{d}{dx}\left[R(x)^{2}\frac{d}{dx}\right]\varphi^{n}_{m}(x)-\varphi^{n}_{m}(x)=-S^{n}_{m}(x),$
(S2)

Using
$R(x)^{2}=R_{0}^{2}\mathrm{e}^{\frac{|x|-L/2}{\delta}}=R_{0}^{2}/(t+t_{0})$,
where $t_{0}=\mathrm{e}^{-\frac{H-L/2}{\delta}}$, the differential operator
can be represented by

$\frac{1}{\kappa^{2}R(x)^{2}}\frac{d}{dx}\left[R(x)^{2}\frac{d}{dx}\right]=\left\\{\begin{array}[]{ll}\displaystyle\frac{(t+t_{0})^{2}}{\kappa^{2}\delta^{2}}\frac{d^{2}}{dt^{2}}&\displaystyle\textrm{for
}|x|>\frac{L}{2},\\\ \\\
\displaystyle\frac{1}{\kappa^{2}}\frac{d^{2}}{dx^{2}}&\displaystyle\textrm{for
}|x|<\frac{L}{2}.\end{array}\right.$ (S3)

Because eq. S2 is an ordinary linear differential equation, the exact solution
can be formally written by

$\varphi^{n}_{m}(x)=\varphi_{\mathrm{p}}(x)+\varphi_{\mathrm{g}}(x),$ (S4)

where $\varphi_{\mathrm{p}}(x)$ is the particular solution of eq. S2, and
$\varphi_{\mathrm{g}}(x)$ is the general solution for the homogeneous
equation. The general solution of the homogeneous differential equation is
given by

$\varphi_{\mathrm{g}}(x)=\left\\{\begin{array}[]{ll}\displaystyle
A_{+}\mathrm{e}^{-\frac{k_{+}(|x|-L/2)}{\delta}}+A_{-}\mathrm{e}^{-\frac{k_{-}(|x|-L/2)}{\delta}}&\displaystyle\textrm{for
}|x|>\frac{L}{2},\\\ \\\ \displaystyle A_{0}\frac{\cosh(\kappa
x)}{\cosh(\kappa L/2)}&\displaystyle\textrm{for
}|x|<\frac{L}{2},\end{array}\right.$ (S5)

where $k_{\pm}=(1\pm\sqrt{1+4\kappa^{2}\delta^{2}})/2$ and $A_{\pm}$ and
$A_{0}$ are adjustable constants. The particular solution
$\varphi_{\mathrm{p}}(x)$ is generally not necessary to satisfies all the
boundary conditions. The boundary conditions are $\varphi^{\prime}(0)=0$,
$\varphi(H)=0$, and the continuity of $\varphi(x)$ and $\varphi^{\prime}(x)$
at $x=L/2$. The adjustable constants $A_{\pm}$ and $A_{0}$ are determined to
satisfies the boundary conditions, and we obtain

$\displaystyle\varphi_{\mathrm{p}}(H)+A_{+}t_{0}^{k_{+}}+A_{-}t_{0}^{k_{-}}$
$\displaystyle=$ $\displaystyle 0,$ (S6)
$\displaystyle\Delta\varphi_{\mathrm{p}}(\frac{L}{2})+A_{+}+A_{-}$
$\displaystyle=$ $\displaystyle A_{0},$ (S7)
$\displaystyle\Delta\varphi_{\mathrm{p}}^{\prime}(\frac{L}{2})-\frac{k_{+}A_{+}+k_{-}A_{-}}{\delta}$
$\displaystyle=$ $\displaystyle A_{0}\kappa\tanh\left(\frac{\kappa
L}{2}\right),$ (S8)

where

$\displaystyle\Delta\varphi_{\mathrm{p}}(\frac{L}{2})=\varphi_{\mathrm{p}}(\frac{L}{2}+0)-\varphi_{\mathrm{p}}(\frac{L}{2}-0),$
(S9)
$\displaystyle\Delta\varphi_{\mathrm{p}}^{\prime}(\frac{L}{2})=\varphi_{\mathrm{p}}^{\prime}(\frac{L}{2}+0)-\varphi_{\mathrm{p}}^{\prime}(\frac{L}{2}-0).$
(S10)

When the salt concentration is high enough for $\kappa L\gg 1$ and
$\kappa\delta\gg 1$, the general solution, $\varphi_{\mathrm{g}}(x)$ is almost
zero except for near the boundary $x=\pm H$ and $\pm L/2$. Therefore, even if
$\psi_{\mathrm{g}}(x)$ is added to satisfy the boundary condition,
$\psi_{\mathrm{g}}(x)$ usually does not significantly affect the charge
distribution, the water flux, and the ionic currents. Direct substitution can
confirm that $\varphi^{1}_{0}(x)=S^{1}_{0}(x)$ and
$\varphi^{0}_{1}(x)=S^{0}_{1}(x)$ are the particular solution of the expanded
Poisson equations. However, these particular solutions do not satisfy the
continuity of ${\varphi^{n}_{m}}^{\prime}(x)$ at $x=\pm L/2$. Thus, we need to
sum $\varphi_{\mathrm{g}}(x)$ to obtain the exact solution, but it is
negligible to calculate the electrokinetic transport coefficients.

For $(n,m)=(2,0)$ and $(1,1)$, the particular solution is

$\varphi_{\mathrm{p}}(x)=\left\\{\begin{array}[]{ll}\displaystyle
S^{n}_{m}(x)+\frac{A^{n}_{m}h\nu}{(\kappa^{2}\delta^{2}-2)}(t+t_{0})^{2}&\displaystyle\textrm{for
}|x|>\frac{L}{2},\\\ \\\ \displaystyle
S^{n}_{m}(x)-\frac{2A^{n}_{m}ht^{\prime}}{\kappa^{2}L^{2}}&\displaystyle\textrm{for
}|x|<\frac{L}{2},\end{array}\right.$ (S11)

where the coefficient $A^{n}_{m}$ is given by

$A^{n}_{m}=\left\\{\begin{array}[]{ll}\displaystyle\frac{\sigma_{0}(2-\sigma_{1})h\nu}{2}&\textrm{for
}(n,m)=(2,0),\\\ \\\ \displaystyle\frac{(2-\sigma_{1})h\nu}{2}&\textrm{for
}(n,m)=(1,1).\end{array}\right.$ (S12)

Because this particular solution does not satisfies the boundary conditions,
we need to add $\varphi_{\mathrm{g}}(x)$ with appropriate $A_{\pm}$ and
$A_{0}$. However, in the case of $\kappa L\gg 1$ and $\kappa\delta\gg 1$, not
only $\varphi_{\mathrm{g}}(x)$ but also $\varphi_{\mathrm{p}}(x)-S^{n}_{m}(x)$
(the second term of eq. S11) are negligible. Thus, we can approximate the
solutions by $\varphi^{2}_{1}(x)=S^{2}_{1}(x)$ and
$\varphi^{1}_{2}(x)=S^{1}_{2}(x)$.

##  VIII Calculation of $\Delta p^{n}_{m}$ and ${J_{i}}^{n}_{m}$ up to $n+m=3$
and conversion to $Q^{n}_{m}$ and $I^{n}_{m}$

Using the definition of $J_{i}$ and the approximated solution of
$\varphi^{n}_{m}(x)$, we obtain ${J_{i}}^{n}_{m}$ up to the third-order as

$\displaystyle{J_{i}}^{0}_{0}$ $\displaystyle=$
$\displaystyle{J_{i}}^{2}_{0}={J_{i}}^{1}_{1}={J_{i}}^{0}_{2}=0,$ (S13)
$\displaystyle{J_{i}}^{1}_{0}$ $\displaystyle=$ $\displaystyle p_{i}h,$ (S14)
$\displaystyle{J_{i}}^{0}_{1}$ $\displaystyle=$ $\displaystyle q_{i}h,$ (S15)
$\displaystyle{J_{i}}^{3}_{0}$ $\displaystyle=$
$\displaystyle\frac{{p_{i}}^{3}}{4}h-3p_{i}h^{2}{f_{i}}^{2}_{0}(H),$ (S16)
$\displaystyle{J_{i}}^{2}_{1}$ $\displaystyle=$
$\displaystyle\frac{q_{i}{p_{i}}^{2}}{4}h-q_{i}h^{2}{f_{i}}^{2}_{0}(H)-2p_{i}h^{2}{f_{i}}^{1}_{1}(H),$
(S17) $\displaystyle{J_{i}}^{1}_{2}$ $\displaystyle=$
$\displaystyle\frac{p_{i}}{4}h-2q_{i}h^{2}{f_{i}}^{1}_{1}(H)-p_{i}h^{2}{f_{i}}^{0}_{2}(H),$
(S18) $\displaystyle{J_{i}}^{0}_{3}$ $\displaystyle=$
$\displaystyle\frac{q_{i}}{4}h-3q_{i}h^{2}{f_{i}}^{0}_{2}(H),$ (S19)

where analytical expressions for ${f_{i}}^{2}_{0}(H)$, ${f_{i}}^{1}_{1}(H)$,
and ${f_{i}}^{0}_{2}(H)$ are as follows:

$\displaystyle{f_{i}}^{2}_{0}(H)$ $\displaystyle=$
$\displaystyle\frac{q_{i}\sigma_{0}(2-\sigma_{1})}{12}h\nu t^{\prime}(1+\nu
t^{\prime})+\frac{h(\alpha_{i}+\sigma_{1}\nu t^{\prime})^{2}}{12}+\frac{\nu
t^{\prime}}{6}(\alpha_{i}+2\nu t^{\prime})\left[\alpha_{i}+2\nu
t^{\prime}+h(\alpha_{i}+\sigma_{1}\nu t^{\prime})\right],$ (S20)
$\displaystyle{f_{i}}^{1}_{1}(H)$ $\displaystyle=$
$\displaystyle\frac{q_{i}h}{12}\left[h\left(\alpha_{i}+\sigma_{1}\nu
t^{\prime}\right)+\nu t^{\prime}\left[(2-\sigma_{1}+6\alpha_{i})(1+\nu
t^{\prime})+2\nu t^{\prime}(3+4\nu t^{\prime})-q_{i}\sigma_{0}h\nu
t^{\prime}(3+2\nu t^{\prime})\right]\right],$ (S21)
$\displaystyle{f_{i}}^{0}_{2}(H)$ $\displaystyle=$ $\displaystyle\frac{1+2\nu
t^{\prime}}{12}.$ (S22)

The pressure difference is calculated as follows:

$\begin{split}\Delta
p&=\frac{2D\eta}{R_{0}b}\mathrm{Pe}+k_{\mathrm{B}}T\rho_{0}\sum_{i,n,m}F(H;\beta_{i}(x){\rho_{i}}^{n}_{m}(x))\frac{(\mathrm{Pe})^{n+1}}{n!}\frac{(\Delta\varphi)^{m}}{m!}-k_{\mathrm{B}}T\rho_{0}\sum_{i,n,m}p_{i}{J_{i}}^{n}_{m}\frac{(\mathrm{Pe})^{n}}{n!}\frac{(\Delta\varphi)^{m}}{m!}.\end{split}$
(S23)

Further calculation yields $\Delta p^{n}_{m}$ up to the third order,

$\displaystyle\Delta p^{0}_{0}$ $\displaystyle=$ $\displaystyle\Delta
p^{2}_{0}=\Delta p^{1}_{1}=\Delta p^{0}_{2}=0,$ (S24) $\displaystyle\Delta
p^{1}_{0}$ $\displaystyle=$
$\displaystyle\sum_{i}p_{i}(1-{J_{i}}^{1}_{0})+\frac{2D\eta}{\rho_{0}k_{\mathrm{B}}TR_{0}b},$
(S25) $\displaystyle\Delta p^{0}_{1}$ $\displaystyle=$
$\displaystyle-\sum_{i}p_{i}{J_{i}}^{0}_{1},$ (S26) $\displaystyle\Delta
p^{3}_{0}$ $\displaystyle=$
$\displaystyle\sum_{i}\left[3F\left(H;\beta_{i}(x){\rho_{i}}^{2}_{0}(x)\right)-p_{i}{J_{i}}^{3}_{0}\right],$
(S27) $\displaystyle\Delta p^{2}_{1}$ $\displaystyle=$
$\displaystyle-\sum_{i}p_{i}{J_{i}}^{2}_{1},$ (S28) $\displaystyle\Delta
p^{1}_{2}$ $\displaystyle=$ $\displaystyle-\sum_{i}p_{i}{J_{i}}^{1}_{2},$
(S29) $\displaystyle\Delta p^{0}_{3}$ $\displaystyle=$
$\displaystyle-\sum_{i}p_{i}{J_{i}}^{0}_{3},$ (S30)

where $F\left(H;\beta_{i}(x){\rho_{i}}^{2}_{0}(x)\right)$ are obtained as

$F\left(H;\beta_{i}(x){\rho_{i}}^{2}_{0}(x)\right)=\frac{(2-\sigma_{1})h\nu
t^{\prime}}{12}\left[\alpha_{i}(\sigma_{1}+6\nu t^{\prime})+4\nu^{2}t^{\prime
2}\right].$ (S31)

Obtained ${J_{i}}^{n}_{m}$ and $\Delta p^{n}_{m}$ are the electrokinetic
coefficients in the expansion of $\mathrm{Pe}$ and $\Delta\varphi$. It is
necessary to convert ${J_{i}}^{n}_{m}$ and $\Delta p^{n}_{m}$ into $Q^{n}_{m}$
and $I^{n}_{m}$ in the expansion of $\Delta P$ and $\Delta\varphi$. $\Delta P$
as a function of $\mathrm{Pe}$ and $\Delta\varphi$ up to the third order is
given by

$\begin{split}\Delta P(\mathrm{Pe},\Delta\varphi)&=\left[\Delta
p^{0}_{1}\Delta\varphi+\Delta
p^{0}_{3}\frac{(\Delta\varphi)^{3}}{6}\right]+\left[\Delta p^{1}_{0}+\Delta
p^{1}_{2}\frac{(\Delta\varphi)^{2}}{2}\right]\mathrm{Pe}+\left[\Delta
p^{2}_{1}\Delta\varphi\right]\frac{(\mathrm{Pe})^{2}}{2}+\left[\Delta
p^{3}_{0}\right]\frac{(\mathrm{Pe})^{3}}{3!}.\end{split}$ (S32)

The Taylor expansion of $\mathrm{Pe}$ as a function of $\Delta P$ is given by

$\mathrm{Pe}=\frac{1}{a_{1}}(\Delta P-a_{0})-\frac{a_{2}}{2{a_{1}}^{3}}(\Delta
P-a_{0})^{2}+\frac{3{a_{2}}^{2}-a_{1}a_{3}}{6{a_{1}}^{5}}(\Delta
P-a_{0})^{3},$ (S33)

where

$\displaystyle a_{0}$ $\displaystyle=$ $\displaystyle\Delta
P|_{\mathrm{Pe}=0}=\Delta p^{0}_{1}\Delta\varphi+\Delta
p^{0}_{3}\frac{(\Delta\varphi)^{3}}{6},$ (S34) $\displaystyle a_{1}$
$\displaystyle=$ $\displaystyle\left.\frac{d\Delta
P}{d\mathrm{Pe}}\right|_{\mathrm{Pe}=0}=\Delta p^{1}_{0}+\Delta
p^{1}_{2}\frac{(\Delta\varphi)^{2}}{2},$ (S35) $\displaystyle a_{2}$
$\displaystyle=$ $\displaystyle\left.\frac{d^{2}\Delta
P}{d\mathrm{Pe}^{2}}\right|_{\mathrm{Pe}=0}=\Delta p^{2}_{1}\Delta\varphi,$
(S36) $\displaystyle a_{3}$ $\displaystyle=$
$\displaystyle\left.\frac{d^{3}\Delta
P}{d\mathrm{Pe}^{3}}\right|_{\mathrm{Pe}=0}=\Delta p^{3}_{0}.$ (S37)

Because

$\mathrm{Pe}=\sum_{n+m=0}^{\infty}Q^{n}_{m}\frac{(\Delta
P)^{n}}{n!}\frac{(\Delta\varphi)^{m}}{m!},$ (S38)

expansion of eq. S33 using eqs. S34–S37 gives

$\displaystyle Q^{1}_{0}$ $\displaystyle=$ $\displaystyle K$ (S39)
$\displaystyle Q^{0}_{1}$ $\displaystyle=$ $\displaystyle-K\Delta p^{0}_{1},$
(S40) $\displaystyle Q^{3}_{0}$ $\displaystyle=$ $\displaystyle-K^{4}\Delta
p^{3}_{0},$ (S41) $\displaystyle Q^{2}_{1}$ $\displaystyle=$
$\displaystyle-K^{3}\Delta p^{2}_{1}+K^{4}\Delta p^{3}_{0}\Delta p^{0}_{1},$
(S42) $\displaystyle Q^{1}_{2}$ $\displaystyle=$ $\displaystyle-K^{2}\Delta
p^{1}_{2}+2K^{3}\Delta p^{2}_{1}\Delta p^{0}_{1}-K^{4}\Delta p^{3}_{0}(\Delta
p^{0}_{1})^{2},$ (S43) $\displaystyle Q^{0}_{3}$ $\displaystyle=$
$\displaystyle-K\Delta p^{0}_{3}+3K^{2}\Delta p^{1}_{2}\Delta
p^{0}_{1}-3K^{3}\Delta p^{2}_{1}(\Delta p^{0}_{1})^{2}+K^{4}\Delta
p^{3}_{0}(\Delta p^{0}_{1})^{3},$ (S44)

where $K=1/\Delta p^{1}_{0}$. Therefore, using eqs. S39–S44, we obtain
$I^{n}_{m}$ given by

$\displaystyle I^{1}_{0}$ $\displaystyle=$
$\displaystyle\sum_{i}q_{i}{J_{i}}^{1}_{0}Q^{1}_{0},$ (S45) $\displaystyle
I^{0}_{1}$ $\displaystyle=$
$\displaystyle\sum_{i}q_{i}\left({J_{i}}^{1}_{0}Q^{0}_{1}+{J_{i}}^{0}_{1}\right),$
(S46) $\displaystyle I^{3}_{0}$ $\displaystyle=$
$\displaystyle\sum_{i}q_{i}\left({J_{i}}^{1}_{0}Q^{3}_{0}+{J_{i}}^{3}_{0}(Q^{1}_{0})^{3}\right),$
(S47) $\displaystyle I^{2}_{1}$ $\displaystyle=$
$\displaystyle\sum_{i}q_{i}\left({J_{i}}^{1}_{0}Q^{2}_{1}+{J_{i}}^{2}_{1}(Q^{1}_{0})^{2}+{J_{i}}^{3}_{0}(Q^{1}_{0})^{2}Q^{0}_{1}\right),$
(S48) $\displaystyle I^{1}_{2}$ $\displaystyle=$
$\displaystyle\sum_{i}q_{i}\left({J_{i}}^{1}_{0}Q^{1}_{2}+{J_{i}}^{1}_{2}Q^{1}_{0}+2{J_{i}}^{2}_{1}Q^{1}_{0}Q^{0}_{1}+{J_{i}}^{3}_{0}Q^{1}_{0}(Q^{0}_{1})^{2}\right),$
(S49) $\displaystyle I^{0}_{3}$ $\displaystyle=$
$\displaystyle\sum_{i}q_{i}\left({J_{i}}^{1}_{0}Q^{0}_{3}+{J_{i}}^{0}_{3}+3{J_{i}}^{1}_{2}Q^{0}_{1}+3{J_{i}}^{2}_{1}(Q^{0}_{1})^{2}+{J_{i}}^{3}_{0}(Q^{0}_{1})^{3}\right).$
(S50)

##  IX Analysis on the previous experimental data

In Ref. 2, the ionic current across a single carbon nanotube is measured. The
current exhibits a quadratic dependence on the applied pressure difference as

$I(\Delta p,\Delta\psi)=\mu\Delta p+\left[G_{0}+G_{2}(\Delta
p)^{2}\right]\Delta\psi,$ (S51)

where $\mu\,\Delta p$ term is negligible compared to the quadratic term
$G_{2}(\Delta p)^{2}\Delta\psi$. The converted electrokinetic coefficients,
$I^{0}_{1}$ and $I^{2}_{1}$, are obtained as follows:

$\displaystyle I^{0}_{1}$ $\displaystyle=$
$\displaystyle\frac{k_{\mathrm{B}}T}{e^{2}\rho_{0}D}\frac{L}{\pi{R_{0}}^{2}}G_{0},$
(S52) $\displaystyle I^{2}_{1}$ $\displaystyle=$
$\displaystyle\frac{2(k_{\mathrm{B}}T)^{3}\rho_{0}}{e^{2}D}\frac{L}{\pi{R_{0}}^{2}}G_{2},$
(S53)

where $\rho_{0}=1\,$M, $T=300\,$K, $D=2.0\times 10^{-9}\,$m2/s, $L=1\,\mu$m,
and $R_{0}=2\,$nm. Table 1 lists the phenomenological parameters, $G_{0}$ and
$G_{2}$, extracted from carbon nanotube 2 and 3 in Supplementary Information
of Ref. 2, and the converted electrokinetic coefficients, $I^{0}_{1}$ and
$I^{2}_{1}$.

Table 1: The electrokinetic parameters in single carbon nanotubes. [2] $\rho_{0}\,$(M) |  $G_{0}$ (nS) |  $G_{2}$ (nS/bar2) |  $I^{0}_{1}$ | $I^{2}_{1}$
---|---|---|---|---
1 |  4.5 |  0.3 |  48 | 4000
1 |  3 |  7.0 |  32 | 93000

## References

  * [1] M. Doi, “Onsager’s variational principle in soft matter,” J. Phys.: Condens. Matter 23, 284118 (2011).
  * [2] A. Marcotte, T. Mouterde, A. Nigués, A. Siria, and L. Bocquet, “Mechanically activated ionic transport across single-digit carbon nanotubes,” Nature Materials 19, 1057-1061 (2020).