Loading [MathJax]/jax/output/SVG/jax.js
Research article

Dynamics of ionic flows via Poisson-Nernst-Planck systems with local hard-sphere potentials: Competition between cations

  • We study a quasi-one-dimensional steady-state Poisson-Nernst-Planck type model for ionic flows through a membrane channel with three ion species, two positively charged with the same valence and one negatively charged. Bikerman's local hard-sphere potential is included in the model to account for ion sizes. The problem is treated as a boundary value problem of a singularly perturbed differential system. Under the framework of a geometric singular perturbation theory, together with specific structures of this concrete model, the existence and uniqueness of solutions to the boundary value problem for small ion sizes is established. Furthermore, treating the ion sizes as small parameters, we derive an approximation of individual fluxes, from which one can further study the qualitative properties of ionic flows and extract concrete information directly related to biological measurements. Of particular interest is the competition between two cations due to the nonlinear interplay between finite ion sizes, diffusion coefficients and boundary conditions, which is closely related to selectivity phenomena of open ion channels with given protein structures. Furthermore, we are able to characterize the distinct effects of the nonlinear interplays between these physical parameters. Numerical simulations are performed to identify some critical potentials which play critical roles in examining properties of ionic flows in our analysis.

    Citation: Peter W. Bates, Jianing Chen, Mingji Zhang. Dynamics of ionic flows via Poisson-Nernst-Planck systems with local hard-sphere potentials: Competition between cations[J]. Mathematical Biosciences and Engineering, 2020, 17(4): 3736-3766. doi: 10.3934/mbe.2020210

    Related Papers:

    [1] Yiwei Wang, Mingji Zhang . Finite ion size effects on I-V relations via Poisson-Nernst-Planck systems with two cations: A case study. Mathematical Biosciences and Engineering, 2024, 21(2): 1899-1916. doi: 10.3934/mbe.2024084
    [2] Hamid Mofidi . New insights into the effects of small permanent charge on ionic flows: A higher order analysis. Mathematical Biosciences and Engineering, 2024, 21(5): 6042-6076. doi: 10.3934/mbe.2024266
    [3] Dawid Czapla, Sander C. Hille, Katarzyna Horbacz, Hanna Wojewódka-Ściążko . Continuous dependence of an invariant measure on the jump rate of a piecewise-deterministic Markov process. Mathematical Biosciences and Engineering, 2020, 17(2): 1059-1073. doi: 10.3934/mbe.2020056
    [4] Natalia L. Komarova . Spatial stochastic models of cancer: Fitness, migration, invasion. Mathematical Biosciences and Engineering, 2013, 10(3): 761-775. doi: 10.3934/mbe.2013.10.761
    [5] Andrew Vlasic, Troy Day . Poisson integral type quarantine in a stochastic SIR system. Mathematical Biosciences and Engineering, 2020, 17(5): 5534-5544. doi: 10.3934/mbe.2020297
    [6] H.Thomas Banks, Shuhua Hu . Nonlinear stochastic Markov processes and modeling uncertainty in populations. Mathematical Biosciences and Engineering, 2012, 9(1): 1-25. doi: 10.3934/mbe.2012.9.1
    [7] Lih-Ing W. Roeger, Z. Feng, Carlos Castillo-Chávez . Modeling TB and HIV co-infections. Mathematical Biosciences and Engineering, 2009, 6(4): 815-837. doi: 10.3934/mbe.2009.6.815
    [8] Emilio N. M. Cirillo, Matteo Colangeli, Adrian Muntean, T. K. Thoa Thieu . A lattice model for active-passive pedestrian dynamics: a quest for drafting effects. Mathematical Biosciences and Engineering, 2020, 17(1): 460-477. doi: 10.3934/mbe.2020025
    [9] Nancy Azer, P. van den Driessche . Competition and Dispersal Delays in Patchy Environments. Mathematical Biosciences and Engineering, 2006, 3(2): 283-296. doi: 10.3934/mbe.2006.3.283
    [10] Pierre Degond, Maxime Herda, Sepideh Mirrahimi . A Fokker-Planck approach to the study of robustness in gene expression. Mathematical Biosciences and Engineering, 2020, 17(6): 6459-6486. doi: 10.3934/mbe.2020338
  • We study a quasi-one-dimensional steady-state Poisson-Nernst-Planck type model for ionic flows through a membrane channel with three ion species, two positively charged with the same valence and one negatively charged. Bikerman's local hard-sphere potential is included in the model to account for ion sizes. The problem is treated as a boundary value problem of a singularly perturbed differential system. Under the framework of a geometric singular perturbation theory, together with specific structures of this concrete model, the existence and uniqueness of solutions to the boundary value problem for small ion sizes is established. Furthermore, treating the ion sizes as small parameters, we derive an approximation of individual fluxes, from which one can further study the qualitative properties of ionic flows and extract concrete information directly related to biological measurements. Of particular interest is the competition between two cations due to the nonlinear interplay between finite ion sizes, diffusion coefficients and boundary conditions, which is closely related to selectivity phenomena of open ion channels with given protein structures. Furthermore, we are able to characterize the distinct effects of the nonlinear interplays between these physical parameters. Numerical simulations are performed to identify some critical potentials which play critical roles in examining properties of ionic flows in our analysis.



    One of the fundamental concerns of physiology is the function of ion channels. Ion channels are approximately cylindrical, hollow proteins with a hole down their middle that provides a controllable path for electrodiffusion of ions (mainly Na+, K+, Ca++ and Cl) through biological membranes, establishing communications among cells and the external environment. This way, ion channels control a wide range of biological functions. The study of ion channels consists of two related major topics: Structures of ion channels and ionic flow properties. With a given structure of an open channel, the main interest is to understand its electrodiffusion property.

    Beyond general electrodiffusion phenomema for electrolytic solutions in bulks or near charged walls, ionic flows through ion channels have more specifics, that is, the study of ionic flows has to take into considerations of global constraints, including boundary conditions in addition to the structure of an ion channel. As already demonstrated by the celebrated works [51,52,53,54,55] of Hodgkin and Huxley for neurons consisting of a population of ion channels and by works in the volume Single-Channel Recording [94] (edited by B. Sakmann and E. Neher) and many works afterwards, the properties of ion channels depend in an extremely rich way on different regions of boundary concentrations and boundary electric potentials. It is the global constraints and the internal structures of membrane channels that make the relevant electrodiffusion properties specific for ion channel problems.

    Taking the structural characteristics into consideration, the basic continuum model for ionic flows is the Poisson-Nernst-Planck (PNP) system, which treats the aqueous medium (within which ions are migrating) as a dielectric continuum [9,10,14,16,19,27,28,29,30,31,40,41,42,43,59,60,68,93]:

    (εr(X)ε0Φ)=e(ns=1zsCs+Q(X)),Jk=0,Jk=1kBTDk(X)Ckμk(X),k=1,2,n, (1.1)

    where XΩ with Ω being a three-dimensional cylindrical-like domain representing the channel, Q(X) is the distribution of the permanent charge along the interior wall of the channel, εr(X) is the relative dielectric coefficient, ε0 is the vacuum permittivity; e is the elementary charge, kB is the Boltzmann constant, T is the absolute temperature; Φ is the electric potential, for the ith ion species, Ci is the concentration, zi is the valence (the number of charges per particle), Ji is the flux density, Di(X) is the diffusion coefficient, and μi(X) is the electrochemical potential.

    Under some reasonable conditions, the PNP system can be derived as a reduced model from molecular dynamics [95], from Boltzmann equation [3], and from variation principles [56,57,58]. The simplest PNP model is the classic PNP (cPNP) system that contains only the ideal component of electrochemical potential. The cPNP system treats ions essentially as point-charges, and neglects ion size effects. It has been simulated [13,14,16,19,40,48] and analyzed [1,4,5,6,32,39,62,69,72,73,85,97,98,102,104,105] to a great extent. However, a major weak point of the cPNP model is that ions are treated as point of charges, which is only reasonable in the extremely dilute setting. Many extremely important properties of ion channels, such as selectivity, rely on ion sizes critically, in particular, for ions that have the same valence, such as sodium Na+ and potassium K+, the main difference is their ionic sizes. Selectivity of ion channels has been a central issue in the study of cells and tissues at least since Hodgkin, Huxley, and Katz discovered the role of Na+ and K+ currents in the action potential [18,49,50]. Selectivity of ion channels has been theorized to occur is in a variety of ways, ranging from structural discussions [23,24,25,26,79,80,108,109] to scaling models [33], kinetic models [45,46] and electrostatic models [20,21,22]. However, all these different theories of selectivity fail to reproduce experimental measurements over a wide range of conditions.

    On the contrary, reduced models, such as the PNP model, using only physical variables, which include the radius of the pore (the channel geometry) and the dielectric constant of the protein, are successful in interpreting the mechanisms that govern selectivity of the L-type Ca++, Na+ channels and also their mutations [9,10,11,12]. It has come as a surprise that a model of selectivity that includes only a few features of the atomic structure has been able to describe the selectivity properties of calcium and sodium channels very well, in all solutions over a wide range of conditions, with only two adjustable parameters using crystal radii of ions. Examining the finite ion size effects on ionic flows, especially for multiple cations with the same valence should provide deep insights into the selectivity phenomenon of ion channels, and this is the motivation of our work.

    Ionic flows through ion channels exhibit extremely rich phenomena, which is why ion channels are nano-scale valves for essentially all activities of living organisms. This is the very reason that it is a great challenge to understand the mechanisms of ion channel functions. To study effects on ionic flows from finite ion sizes, one has to consider excess (beyond the ideal) component in the electrochemical potential. One way is to include uncharged hard-sphere potentials to partially account for ion size effects. Physically, this means that each ion is approximated as a hard-sphere with its charges at the center of the sphere. Both local and nonlocal models for hard-sphere potentials were introduced for this purpose [8,91,92]. Nonlocal models give the hard-sphere potentials as functionals of ion concentrations while local models depend pointwise on ion concentrations. The PNP models with ion sizes have been investigated analytically [2,7,61,63,71,76,77,99,106] for two ion species, one positively charged and one negatively charged and computationally [41,42,43,47,56,57,58,68,84,88,96,107] for ion channels and have shown great success. Existence and uniqueness of minimizers and saddle points of the free-energy equilibrium formulation with ionic interaction have also been mathematically analyzed [34,70].

    For mathematical analysis, the challenge lies in the fact that specific dynamics depend on complicated nonlinear interplays of multiple physical parameters such as boundary concentrations and potentials, diffusion coefficients, ion sizes and valences, permanent charge distributions, etc. In general, there is no hope to have explicit solution formulae for such a complicated problem even with simple boundary values. The recent development in analyzing classic PNP models [32,72,73] sheds some lights on the voltage-current relationship in simplified settings. This development is based heavily on modern invariant manifold theory of nonlinear dynamical systems, particularly, the geometric theory of singular perturbations. Together with specific structures to PNP models, we are able to, far beyond the existence results, obtain a more or less explicit approximation formula for solutions from which one can further study the qualitative properties of ionic flows and extract concrete information directly related to biological measurements.

    In this work, we will study the PNP model with three ion species, two positively charged with the same valence and one negatively charged. Bikerman's local hard-sphere model is included to account for finite ion size effects. Of particular interest is the competition between two cations due to finite ion sizes, which provides useful and deep insights into the selectivity phenomena for open ion channels with given protein structures. Our analysis is based on a further reduction of PNP models. On the basis that ion channels have narrow cross sections relative to their lengths, PNP systems defined on three-dimensional ion channels are further reduced to quasi-one-dimensional models first proposed in [82] and, for a special case, the reduction is rigorously justified in [75]. A quasi-one-dimensional (time-evolution) PNP type model for ion flows of n ion species is

    1A(X)ddX(εr(X)ε0A(X)dΦdX)=e(nj=1zjCj(X)+Q(X)),Cit+dJidX=0,Ji=1kBTDi(X)A(X)Ci(X)dμidX,i=1,2,,n, (1.2)

    where A(X) is the area of the cross-section of the channel at location X. The boundary conditions are, for i=1,2,,n,

    Φ(t,0)=V,  ci(t,0)=Li>0;Φ(t,1)=0,  ci(t,1)=Ri>0. (1.3)

    For ion channels, an important characteristic is the I-V relation. Given a solution of the boundary value problem (1.2)–(1.3), the current is

    I=nk=1zkJk, (1.4)

    where zkJk is the individual flux of charge of the kth ion species. For fixed boundary concentrations Lk's and Rk's, Jk's depend on V only and formula (1.4) provides a dependence of the current I(V;ε) on the voltage V. In terms of applications, it is important to study properties of individual fluxes Jk due to the fact that most experiments (with some exceptions) can only measure the total current I while individual fluxes contain much more information on channel functions. This is another reason that we mainly focus on the study of individual fluxes in this work. We would like to point out that in general, the I-V relation is not unique for the PNP system with nonzero permanent charge Q(x), even for the classic PNP system with two ion species [32]. However, with the permanent charge Q(x) being zero, the I-V relation is unique, that is, the cPNP system has a unique solution [72]. In our work, we assume the permanent charge to be zero, and view the volume of the ion species as a regular perturbation to the system, and since the linearized system at zero ion-size is nondegenerate, the solution to the PNP system is unique; but the dynamics of the ionic flows is still rich due to the nonlinear interplay between distinct physical parameters.

    We comment that our results, for the relatively simple setting and assumptions of our model, are rigorous. We believe these results will provide useful insights for numerical and even experimental studies of ionic flows through membrane channels. At the same time, we point out that the quasi-one-dimensional PNP model and the local hard-sphere model adopted in this paper are rather simple. Aside the fact that they will miss the three-dimensional features of the problem, a major weakness is the omission of the excess electrostatic component in the excess potentials. As a result, important phenomena such as charge inversion and layering may not be detected by this model.

    The rest of this paper is organized as follows. In section 2, we set up our problem with further assumptions. In section 3, the existence and (local) uniqueness result for the boundary value problem is established under the framework of the geometric singular perturbation theory. Based on the analysis, treating the ion sizes as small parameters, approximations of individual fluxes are derived, from which the effects on ionic flows from ion sizes are analyzed in details. This leads to our main focus, the competition between two cations due to the nonlinear interplay between finite ion sizes, diffusion coefficients and boundary conditions, studied in section 4. Some concluding remarks are stated in section 5. A detailed proof of Lemma 3.5 is provided in Appendix A.

    The electrochemical potential μi(X) for the ith ion species consists of the ideal component μidi(X) and the excess component μexi(X): μi(X)=μidi(X)+μexi(X), where

    μidi(X)=zieΦ(X)+kBTlnCi(X)C0 (2.1)

    with some characteristic number density C0. The classical PNP system takes into consideration of the ideal component μidi(X) only. This component reflects the collision between ion particles and the water molecules.

    It has been accepted that the classic PNP system is a reasonable model in, for example, the dilute case under which the ion particles can be treated as point particles and the ion-to-ion interaction can be more or less ignored.

    The excess chemical potential μexi(X) accounts for finite sizes of ions, which is critically important for selectivitywhich ion a channel prefersof ion channel. While channel protein structures play an important role for channel selectivity, the electrodiffusion of ions is clearly another critical component. For different ions with the same valence, their sizes are a crucial factor for the selectivity. The excess potential μexi(X) consists of two components: the (uncharged) hard-sphere component μHSi(X) and the excess electrostatic component μESi. While there are very successful models for μHSi(X), modeling of μESi(X) is itself extremely challenging and is still a very active research area in liquid-state theory of chemistry and physics [34,41,90,91].

    In this paper, we will take the following Bikerman's local hard-sphere model to approximate μexi(X)

    1kBTμbiki(X)=ln(1nj=1νjCj(X))=νnj=1λjCj(X)+o(ν), (2.2)

    where νj is the volume of the jth ion species. In particular, in the following analysis, we will take νn=ν so that λn=1.

    The main goal of this paper is to examine the qualitative effect of ion sizes via the steady-state boundary value problem of (1.2)–(1.3) with the local hard-sphere model (2.2) for the excess potential. For definiteness, we will take the following settings:

    (A1) We consider three ion species (n=3) with z1=z2=z3=1.

    (A2) The permanent charge is set to be zero: Q(X)=0.

    (A3) For the electrochemical potential μi, in addition to the ideal component μidi, we also include the local hard-sphere potential μbiki in (2.2).

    (A4) The relative dielectric coefficient and the diffusion coefficient are constants, that is, εr(X)=εr and Di(X)=Di.

    Remark 2.1. In the study of ion channel problems, the selectivity of cations Na+ and K+ is extremely important. This is exactly the reason that we choose z1=z2=1 together with z3=1 for the anion Cl.

    In the sequel, we will assume (A1)–(A4). Under the assumptions (A1)–(A4), the steady-state system of (1.2) reads

    1A(X)ddX(εrε0A(X)dΦdX)=e3j=1zjCj(X),dJidX=0,Ji=1kBTDiA(X)Ci(X)dμidX,i=1,2,3. (2.3)

    The boundary conditions become, for i=1,2,3,

    Φ(0)=V,  ci(0)=Li>0;Φ(1)=0,  ci(1)=Ri>0. (2.4)

    We first make a dimensionless rescaling following [39]. Set C0=max{Li,Ri:i=1,2} and let

    ε2=εrε0kBTe2l2C0, x=Xl, h(x)=A(X)l2, Di=lC0Di; ϕ(x)=ekBTΦ(X),ci(x)=Ci(X)C0, Ji=JiDi; V=ekBTV, Li=LiC0; Ri=RiC0. (2.5)

    We would like to point out that the dimensionless parameter ε defined in (2.5) as ε=1lεrε0kBTe2C0 is directly related to the ratio κD/l, where κD=εrε0kBTj(zje)2Cj is the Debye length; in particular, ε=κD/l when z2j=1 and Cj=C0. Typically, the parameter ε is small due to the fact that the two variables l, the length of the channel, and C0, some characteristic number density could be very large (for many cases, the value of ε is of order O(103)).

    The boundary value problem (2.3)–(2.4) becomes

    ε2h(x)ddx(h(x)ddxϕ)=(c1+c2c3),dc1dx+c1dϕdx+c1(x)kBTddxμbik1(x)=J1h(x),dc2dx+c2dϕdx+c2(x)kBTddxμbik2(x)=J2h(x),dc3dxc3dϕdx+c3(x)kBTddxμbik3(x)=J3h(x),dJkdx=0, k=1,2,3 (2.6)

    with the boundary condition

    ϕ(0)=V,  ck(0)=Lk>0;  ϕ(1)=0,  ck(1)=Rk>0. (2.7)

    Remark 2.2. We will take h(x)=1 over the whole interval [0,1] in our analysis. This is because for ion channels with zero permanent charge, it turns out that the variable h(x) contributes through an average, explicitly through the factor 110h1(x)dx (for example, the individual flux will be DkJk10h1(x)dx, see [71]), which does not affect our analysis of qualitative properties of ionic flows.

    Denote

    C=(c1c2c3),J=(J1J2J3),Z=diag(1,1,1),N=(c1cμbik1c2cμbik2c3cμbik3)=(c1c1μbik1c1c2μbik1c1c3μbik1c2c1μbik2c2c2μbik2c2c3μbik2c3c1μbik3c3c2μbik3c3c3μbik3).

    Then, system (2.6) can be rewritten as

    ε2d2ϕdx2=(c1+c2c3),ddxC=dϕdx(I+N)1ZC(x)(I+N)1J,ddxJ=0. (2.8)

    Recall that the boundary condition is now

    ϕ(0)=V,  ck(0)=Lk;  ϕ(1)=0,  ck(1)=Rk. (2.9)

    We will rewrite system (2.8) into a standard form for singularly perturbed systems and convert the boundary value problem (2.8) and (2.9) to a connecting problem.

    Denote the derivative with respect to x by overdot and introduce u=ε˙ϕ and τ=x. System (2.8) becomes

    ε˙ϕ=u,ε˙u=(c1+c2c3),ε˙C=u(I+N)1ZCε(I+N)1J,˙J=0,˙τ=1. (3.1)

    System (3.1) will be treated as a dynamical system of phase space R9 with state variables (ϕ,u,C,J,τ).

    For ε>0, the rescaling x=εξ of the independent variable x gives rise to the fast system,

    ϕ=u,u=(c1+c2c3),C=u(I+N)1ZCε(I+N)1J,J=0,τ=ε, (3.2)

    where prime denotes the derivative with respect to ξ.

    Let BL and BR be the subsets of the phase space R9 defined by

    BL={(V,u,L1,L2,L3,J1,J2,J3,0)R9:arbitrary u,J1,J2,J3},BR={(0,u,R1,R2,R3,J1,J2,J3,1)R9:arbitrary u,J1,J2,J3}. (3.3)

    Then the original boundary value problem is equivalent to a connecting problem, namely, finding a solution of (3.1) or (3.2) from BL to BR.

    In what follows, instead of studying the boundary value problem, we will consider the equivalent connecting problem for system (3.1) or (3.2) and construct its solution from BL to BR. The construction process involves two main steps: the first step is to construct a singular orbit to the connecting problem, and the second step is to apply geometric singular perturbation theory to show that there is a unique solution near the singular orbit for small ε>0.

    Following the idea in [32,71,72,73], we first construct a singular orbit on [0,1] that connects BL to BR. In our coming analysis, we assume the so-called electroneutrality conditions

    L1+L2=L3 and R1+R2=R3. (3.4)

    By setting ε=0 in (3.1), we obtain the so-called slow manifold

    Z={u=0,c1+c2c3=0}. (3.5)

    By setting ε=0 in (3.2), we get the limiting fast system

    ϕ=u,u=(c1+c2c3),C=u(I+N)1ZC,J=0,τ=0. (3.6)

    Note that the slow manifold Z is precisely the set of equilibria of (3.6). It follows from [32,61,63,71,72,73,76] that

    Lemma 3.1. For system (3.6), the slow manifold Z is normally hyperbolic.

    We would like to point out that the normal hyperbolicity of the slow manifold Z is crucial in our analysis, which guarantees the existence of the stable manifold and the unstable manifold that will be discussed below (see [32,35] for more details for the theory of normally hyperbolic invariant manifolds).

    For ν>0 small, treating (3.6) as a regular perturbation of that with ν=0, we look for solutions Γ(ξ;ν)=(ϕ(ξ;ν),u(ξ;ν),c1(ξ;ν),c2(ξ;ν),c3(ξ;ν),J1(ν),J2(ν),J3(ν),τ) of system (3.6) of the form

    ϕ(ξ;ν)=ϕ0(ξ)+ϕ1(ξ)ν+o(ν),u(ξ;ν)=u0(ξ)+u1(ξ)ν+o(ν),ck(ξ;ν)=ck0(ξ)+ck1(ξ)ν+o(ν),Jk(ν)=Jk0+Jk1ν+o(ν). (3.7)

    Lemma 3.2. For ν>0 small, one has, up to the first order in ν,

    (I+N)1=(100010001)(λ1c10λ2c10c10λ1c20λ2c20c20λ1c30λ2c30c30)ν+o(ν). (3.8)

    Proof. A direct computation leads to the result. We omit it here.

    Substituting (3.7) into system (3.6), together with (3.8), by careful calculations, we obtain the zeroth order limiting fast system in ν:

    ϕ0=u0,u0=(c10+c20c30),c10=u0c10,c20=u0c20,c30=u0c30,J10=J20=J30=0,τ=0, (3.9)

    and the first order limiting fast system in ν:

    ϕ1=u1,u1=(c11+c21c31),(c11c21c31)=(u0c11+u1c10u0c21+u1c20u0c31u1c30)+(λ1c10λ2c10c10λ1c20λ2c20c20λ1c30λ2c30c30)(u0c10u0c20u0c30),J11=J21=J31=0,τ=0. (3.10)

    Let ML be the collection of orbits from BL in forward time under the flow of system (3.6) and MR be the collection of orbits from BR in backward time under the flow of system (3.6). Let Ws(Z) be the stable manifold of Z that consists of points approaching Z in forward time, and let Wu(Z) be the unstable manifold of Z that consists of points approaching Z in backward time. Then, for a singular orbit connecting BL to BR, the boundary layer at x=0 must lie in NL=MLWs(Z) and the boundary layer at x=1 must lie in NR=MRWu(Z). Let ω(NL) be the ω-limit set of NL, and α(NR) be the α-limit set of NR, which are defined by

    ω(NL)={ϕ0()+ϕ1()ν+o(ν),0,c10()+c11()ν+o(ν),c20()+c21()ν+o(ν),c30()+c31()ν+o(ν),J1(ν),J2(ν),J3(ν),0},α(NR)={ϕ0()+ϕ1()ν+o(ν),0,c10()+c11()ν+o(ν),c20()+c21()ν+o(ν),c30()+c31()ν+o(ν),J1(ν),J2(ν),J3(ν),0}.

    Proposition 3.3. Assume that ν0 is small. Under electronutrality conditions (3.4), one has

    (ⅰ) The stable manifold Ws(Z) intersects BL transversally at points

    (V,ul0+ul1ν+o(ν),L1,L2,L3,J1(ν),J2(ν),J3(ν),0)

    and the ω-limit set of NL=MLWs(Z) is

    ω(NL)={ϕL0+ϕL1ν+o(ν),0,cL10+cL11ν+o(ν),cL20+cL21ν+o(ν),cL30+cL31ν+o(ν),J1(ν),J2(ν),J3(ν),0}

    where Ji(ν)=Ji0+Ji1ν+o(ν),i=1,2,3, can be arbitrary, and

    ϕL0=VcL10=L1,cL20=L2,cL30=L3,ul0=2L3,ϕL1=cL11=cL21=cL31=uL1=0,

    (ⅱ) The unstable manifold Wu(Z) intersects BR transversally at points

    (0,ur0+ur1ν+o(ν),R1,R2,R3,J1(ν),J2(ν),J3(ν),1)

    and the α-limit set of NR=MRWu(Z) is

    α(NR)={(ϕR0+ϕR1ν+o(ν),0,cR10+cR11ν+o(ν),cR20+cR21ν,cR30+cR31ν+o(ν),J1(ν),J2(ν),J3(ν),1)}

    where Ji(ν)=Ji0+Ji1ν+o(ν), i=1,2,3, can be arbitrary, and

    ϕR0=0,cR10=R1,cR20=R2,cR30=R3,ur0=2R3,ϕR1=cR11=cR21=ur1=0.

    Proof. The result for system (3.9), the zeroth order system in ν, has been obtained in [6], we will not repeat it here. For system (3.10), the first order system in ν, it has the following four nontrivial first integrals:

    H1=ϕ1+c11c10+λ1c10+λ2c20+c30,H2=ϕ1+c21c20+λ1c10+λ2c20+c30,H3=ϕ1+c31c30+λ1c10+λ2c20+c30,H4=u0u1c11c21c3112(λ1c210+λ2c220+c230)λ1+λ22c10c20+(λ11)c10c30ϕ0+(λ21)c20c30ϕ0.

    We now establish the results for ϕL1,cL11,cL21,cL31 and ul1 for system (3.10). Those for ϕR1,cR11,cR21 and cR31 can be established in a similar way.

    Note that ϕ(0)=c11(0)=c21(0)=c31(0)=0. Using the integrals H1,H2 and H3, one has

    ϕ1+c11c10+λ1c10+λ2c20+c30=λ1L1+λ2L2+L3,ϕ1+c21c20+λ1c10+λ2c20+c30=λ1L1+λ2L2+L3,ϕ1+c31c30+λ1c10+λ2c20+c30=λ1L1+λ2L2+L3.

    Let ξ, careful calculations give

    ϕL1=0,cL11=cL10(λ1L1+λ2L2+L3λ1cL10λ2cL20cL30ϕL1),cL21=cL20(λ1L1+λ2L2+L3λ1cL10λ2cL20cL30ϕL1),cL31=cL30(λ1L1+λ2L2+L3λ1cL10λ2cL20cL30+ϕL1).

    Under electroneutrality conditions (3.4), one has cL10=L1, cL20=L2 and cL30=L3. It follows that cL11=cL21=cL31=0. In view of H4(0)=H4(), together with the above analysis, one has ul1=0. This completes the proof.

    We now construct the regular layer Λ on Z that connects ω(NL) and α(NR). Following the ideas in [32,71,72,73], we make a rescaling u=εp and c1+c2c3=εq in system (3.1). In term of the new variables, system (3.1) becomes

    ˙ϕ=p,ε˙p=q,ε˙q=p(1,1,1)(I+N)1ZC+(1,1,1)(I+N)1J,˙C=p(I+N)1ZC(I+N)1J,˙J=0,˙τ=1, (3.11)

    where c3=c1+c2+εq.

    It is again a singular perturbation problem and its limiting slow system is

    ˙ϕ=p,q=0,0=p(1,1,1)(I+N)1ZC+(1,1,1)(I+N)1J,˙C=p(I+N)1ZC(I+N)1J,˙J=0,˙τ=1, (3.12)

    where c3=c1+c2.

    For system (3.12), the slow manifold is

    S={p=(1,1,1)(I+N)1J(1,1,1)(I+N)1ZC,q=0}.

    Therefore, the limiting slow system on S is

    ˙ϕ=(1,1,1)(I+N)1J(1,1,1)(I+N)1ZC,˙C=(1,1,1)(I+N)1J(1,1,1)(I+N)1ZC(I+N)1ZC(I+N)1J,˙J=0,˙τ=1, (3.13)

    where c3=c1+c2.

    As for the regular layer problem, we look for solutions of (3.13) of the form

    ϕ(x)=ϕ0(x)+ϕ1(x)ν+o(ν),ck(x)=ck0(x)+ck1(x)ν+o(ν),Jk=Jk0+Jk1ν+o(ν) (3.14)

    to connect ω(NL) and α(NL) given in Proposition 3.3; in particular, for j=0,1,

    (ϕj(0),c1j(0),c2j(0))=(ϕLj,cL1j,cL2j),(ϕj(1),c1j(1),c2j(1))=(ϕRj,cR1j,cR2j).

    For convenience, we define quantities Tc0, Tc1, Tm0, Tm1, F1, F2 and a function A(x) as, for k=0,1,

    Tck=J1k+J2kJ3k,Tmk=J1k+J2k+J3k,F1=J20L1J10L2J10+J20,F2=λ1J10+λ2J20J10+J20,A(x)=1Tm02L3x. (3.15)

    Substituting (3.14) into system (3.13), together with (3.8), by careful calculations, we obtain the zeroth order limiting slow system in ν

    ˙ϕ0=Tc02(c10+c20),˙c10=Tc02(c10+c20)c10J10,˙c20=Tc02(c10+c20)c20J20,˙J10=˙J20=˙J30=0,˙τ=1, (3.16)

    and the first order limiting slow system in ν

    ˙ϕ1=Tc0(c11+c21)2(c10+c20)2Tc12(c10+c20),˙c11=Tc0(c11c20c21c10)2(c10+c20)2+((λ1+1)Tm0+2(λ2λ1)J20)c2102(c10+c20)+(λ2+1)Tm0+2(λ1λ2)J102(c10+c20)c10c20+Tc1c102(c10+c20)J11,˙c21=Tc0(c11c20c21c10)2(c10+c20)2+((λ2+1)Tm0+2(λ1λ2)J10)c2202(c10+c20)+(λ1+1)Tm0+2(λ2λ1)J202(c10+c20)c10c20+Tc1c202(c10+c20)J21. (3.17)

    For the zeroth order system (3.16) (see [6] for details), one has

    Lemma 3.4. Under electroneutrality conditions (3.4), there is a unique solution (ϕ0(x),c10(x),c20(x),J10,J20,J30,τ(x)) of (3.16) such that (ϕ0(0),c10(0),c20(0),τ(0))=(ϕL0,cL10,cL20,0,) and (ϕ0(1),c10(1),c20(1),τ(1))=(ϕR0,cR10,cR20,1), where ϕL0,cL10,cL20,ϕR0,cR10, and cR20 are given in Proposition 3.3. It is given by

    ϕ0(x)=V+Tc0Tm0lnA(x),c10(x)=F1ATc0Tm0(x)+J10L3J10+J20A(x),c20(x)=F1ATc0Tm0(x)+J20L3J10+J20A(x),τ(x)=x,J10=f0(L3,R3)f1(L3,R3;V)(L1R1eV),J20=f0(L3,R3)f1(L3,R3;V)(L2R2eV),J30=f0(L3,R3)(lnL3ln(R3eV)), (3.18)

    where

    f0(L3,R3)=L3R3lnL3lnR3,f1(L3,R3;V)=lnL3ln(R3eV)L3R3eV. (3.19)

    Note that, for x(0,1), under (3.4), one has A(x)=1x+R3L3x>0. For convenience, we define functions Mk(x)=M1(L1,L2,R1,R2;λ1,λ2;x) as, k=1,2,3,4,

    M1(x)=[F1Tc0Tm0(A1(x)1)J10L3J10+J20(ATc0Tm0(x)1)](Tm1Tm0(λ2λ1)F1Tc02J30L3(1+F2)2),M2(x)=[F1Tc0Tm0(A1(x)1)+J20L3J10+J20(ATc0Tm0(x)1)]((λ2λ1)F1Tc02J30Tm1Tm0+L3(1+F2)2),M3(x)=F1L3Tm0[(λ2λ1)Tc0J10J10+J20(1+Tc02J30)+(1+F2)(Tm0+Tc02)](A(x)1),M4(x)=F1L3Tm0[(λ2λ1)Tc0J20J10+J20(1Tc02J30)+(1+F2)(Tm0+Tc02)](A(x)1). (3.20)

    For the first order system (3.17), we have (see Appendix A for a detailed proof)

    Lemma 3.5. Under electroneutrality conditions (3.4), there is a unique solution (ϕ1(x),c11(x),c21(x),J11,J21,J31,τ(x)) of (3.17) such that (ϕ1(0),c11(0),c21(0),τ(0))=(ϕL1,cL11,cL21,0)=(0,0,0,0), and (ϕ1(1),c11(1),c21(1),τ(1))=(ϕR1,cR11,cR21,1)=(0,0,0,1), where ϕL1, cL11, cL21, ϕR1, cR11 and cR21 are given in Proposition 3.3. It is given by

    ϕ1(x)=Tc0Tm0(L3(1+F2)2+(λ2λ1)F1Tc02J30Tm1Tm0)(A1(x)1)+L3(1+F2)Tc02Tm0(A(x)1)(λ2λ1)F1Tm02J30(ATc0Tm0(x)1)+Tc1Tm0Tm1Tc0(Tm0)2lnA(x),c11(x)=ATc0Tm0(x){F1(Tc0Tm1Tc1Tm0)(Tm0)2lnA(x)(1+F2)J10L232(J10+J20)(A2+Tc0Tm0(x)1)+M1(x)M3(x)+L3J10+J20(J10(Tm1Tc0Tc1Tm0)2Tm0(J10+J20)+J11)(A1+Tc0Tm0(x)1)+(λ2λ1)F21(1+Tc02J30)(ATc0Tm0(x)1)},c21(x)=ATc0Tm0(x){F1(Tc0Tm1Tc1Tm0)(Tm0)2lnA(x)(1+F2)J20L232(J10+J20)(A2+Tc0Tm0(x)1)+M2(x)+M4(x)+L3J10+J20(J20(Tm1Tc0Tc1Tm0)2Tm0(J10+J20)+J21)(AI0+T0T0(x)1)(λ2λ1)F21(1+Tc02J30)(ATc0Tm0(x)1)},J11=[L3J10+J20(A1+Tc0Tm0(1)1)]1{(1+F2)J10L232(J10+J20)(A2+Tc0Tm0(1)1)+F1(Tc0Tm1Tc1Tm0)(Tm0)2lnA(1)+J10(Tm1Tc0Tc1Tm0)L32Tm0(J10+J20)2(A1+Tc0Tm0(1)1)+(λ2λ1)F21(1+Tc02J30)(ATc0Tm0(1)1)+M1(1)M3(1)},J21=[L3J10+J20(A1+Tc0Tm0(1)1)]1{(1+F2)J20L232(J10+J20)(A2+Tc0Tm0(1)1)+F1(Tc0Tm1Tc1Tm0)(Tm0)2lnA(1)+J20(Tm1Tc0Tc1Tm0)L32Tm0(J10+J20)2(A1+Tc0Tm0(1)1)(λ2λ1)F21(1+Tc02J30)(ATc0Tm0(1)1)+M2(1)+M4(1)},J31=Tm1Tc12,

    where J10,J20,J30,Tm0,Tc0 and Mk(x), k=1,2,3,4 are given in (3.15), (3.18), (3.20), and

    Tm1=L3[(λ2λ1)F1Tc0J30(A1Tc0Tm0(1)1)+L3(1+F2)(A2(1)1)],Tc1=Tm0lnA(1)[Tc0Tm0(Tm1Tm0L3(1+F2)2(λ2λ1)F1Tc02J30)(A1(1)1)L3(1+F2)Tc02Tm0(A(1)1)+(λ2λ1)F1Tc02J30(ATc0Tm0(1)1)]+Tm1Tc0Tm0.

    The next result can be justified directly.

    Proposition 3.6. Assume electroneutrality conditions (3.4). Fixing V and treating both the zeroth order approximations Jk0's (in ν) and first order approximations Jk1's (in ν) as functions of boundary concentrations Li and Ri, one has Jk0 is homogeneous of degree one and Jk1 is homogeneous of degree two, that is, for any α>0,

    Jk0(αL1,αL2,αL3;αR1,αR2,αR3)=αJk0(L1,L2,L3;R1,R2,R3),
    Jk1(αL1,αL2,αL3;αR1,αR2,αR3)=α2Jk1(L1,L2,L3;R1,R2,R3).

    Remark 3.7. This interesting observation actually provides a nice way to adjust the effects on ionic flows (in particular, for the individual flux) from finite ion size by suitably controlling its boundary concentrations. Furthermore, in our later study of the function P(V)=D1J11D2J21, which provides detailed information that is related to the selectivity of cations due to electrodiffusion properties of ions with a given protein structures, this scaling law can either reduce or enhance the preference of the ion channel for cations.

    The slow orbit, up to O(ν),

    Λ(x;ν)=(ϕ0(x)+ϕ1(x)ν,c10(x)+c11(x)ν,c20(x)+c21(x)ν,J10+J11ν,J20+J21ν,J30+J31ν,τ(x)) (3.21)

    given in Lemmas 3.4 and 3.5 connects ω(NL) and α(NR). Let ˉML (resp., ˉMR) be the forward (resp., backward) image of ω(NL) (resp., α(NR)) under the slow flow (3.13) on the seven-dimensional slow manifold S. Following the idea in [71,72], we have

    Proposition 3.8. There exists ν0>0 small depending on boundary conditions so that, if 0νν0, then, on the seven-dimensional slow manifold S, ˉML and ˉMR intersects transversally along the unique orbit Λ(x;ν) given in (3.21).

    Figure 1.  A singular orbit Γ0ΛΓ1 on [0,1]: a boundary layer Γ0 at τ=0, a regular layer Λ on Z from τ=0 to τ=1, and a boundary layer Γ1 at τ=1.

    We have constructed a unique singular orbit on [0, 1] that connects BL to BR. It consists of two boundary layer orbits Γ0 from the point

    (V,ul0+ul1ν+o(ν),L1,L2,L3,J10+J11ν+o(ν),J20+J21ν+o(ν),J30+J31ν+o(ν),0)BL

    to the point (ϕL,0,cL1,cL2,cL3,J1,J2,J3,0)ω(NL)Z and Γ1 from the point (ϕR,0,cR1,cR2,cR3,J1,J2,J3,1)z1(NR)Z to the point

    (0,ur0+ur1ν+o(ν),R1,R2,R3,J1,J2,J3,1)BR,

    and a regular layer Λ on Z that connects the two landing points

    (ϕL,0,cL1,cL2,cL3,J1,J2,J3,0)ω(NL)  and  (ϕR,0,cR1,cR2,cR3,J1,J2,J3,1)α(NR)

    of the two boundary layers.

    We now establish the existence of a solution of (2.8)–(2.9) near the singular orbit constructed above which is a union of two boundary layers and one regular layer Γ0ΛΓ1. The proof follows the same line as that in [6,32,71,72,73] and the main tool used is the Exchange Lemma [64,65,66,100] of the geometric singular perturbation theory.

    Theorem 3.9. Let Γ0ΛΓ1 be the singular orbit of the connecting problem system (3.1) associated to BL and BR in system (3.3). Let ν0>0 be as in Proposition 3.8. Then, there exists ε0>0 small (depending on the boundary conditions and ν0) so that, if 0νν0 and 0<εε0, then the boundary value problem (2.8)–(2.9) has a unique smooth solution near the singular orbit Γ0ΛΓ1.

    Proof. Let ν0 be as in Proposition 3.8. For 0νν0, we define ul=ul0+ul1ν, J1(ν)=J10+J11ν, J2(ν)=J20+J21ν and J3(ν)=J30+J31ν. Fix δ>0 small to be determined. Let

    BL(δ)={(V,u,L1,L2,L3,J1,J2,J3,0)R9:|uul|<δ,|JiJi(ν)|<δ}.

    For ε>0, let ML(ε,δ) be the forward trace of BL(δ) under the flow of system (3.1) or equivalently of system (3.2) and let MR(ε) be the backward trace of BR. To prove the existence and uniqueness statement, it suffices to show that ML(ε,δ) intersects MR(ε) transversally in a neighborhood of the singular orbit Γ0ΛΓ1. The latter will be established by an application of the Exchange Lemma.

    Notice that dimBL(δ) = 4. It is clear that the vector field of the fast system (3.2) is not tangent to BL(δ) for ε0, and hence, dimML(ε,δ) = 5. We next apply Exchange Lemma to track ML(ε,δ) in the vicinity of Γ0ΛΓ1. First of all, the transversality of the intersection BL(δ)Ws(Z) along Γ0 in Proposition 3.3 implies the transversality of intersection ML(0,δ)Ws(Z). Secondly, we have also established that dimω(NL)=dimNL1=3 in Proposition 3.3 and that the limiting slow flow is not tangent to ω(NL) in Section 3.2. Under these conditions, the Exchange Lemma [64,65,66,100] states that there exist ρ>0 and ε1>0 so that, if 0<εε1, then ML(ε,δ) will first follow Γ0 toward ω(NL)Z, then follow the trace of ω(NL) in the vicinity of Λ toward {τ=1}, leave the vicinity of Z, and, upon exit, a portion of ML(ε,δ) is C1O(ε)-close to Wu(ω(NL)×(1ρ,1+ρ)) in the vicinity of Γ1. Note that dimWu(ω(NL)×(1ρ,1+ρ))=dimML(ε,δ)=5.

    It remains to show that Wu(ω(NL)×(1ρ,1+ρ)) intersects MR(ε) transversally since ML(ε,δ) is C1O(ε)-close to Wu(ω(NL)×(1ρ,1+ρ)). Recall that, for ε=0, MR intersects Wu(Z) transversally along NR (Proposition 3.3); in particular, at γ1:=α(Γ1)α(NR)Z, we have

    Tγ1MR=Tγ1α(NR)Tγ1Wu(γ1)span{Vs}

    where, Tγ1Wu(γ1) is the tangent space of the one-dimensional unstable fiber Wu(γ1) at γ1 and the vector VsTγ1Wu(Z) (the latter follows from the transversality of the intersection of MR and Wu(Z)). Also,

    Tγ1Wu(ω(NL)×(1ρ,1+ρ))=Tγ1(ω(NL)1)span{Vτ}Tγ1Wu(γ1)

    where the vector Vτ is the tangent vector to the τ-axis as the result of the interval factor (1ρ,1+ρ). From Proposition 3.8, ω(NL)1 and α(NR) are transversal on Z{τ=1}. Therefore, at γ1, the tangent spaces Tγ1MR and Tγ1Wu(ω(NL)×(1ρ,1+ρ)) contain seven linearly independent vectors: Vs, Vτ, Tγ1Wu(γ1) and the other four from Tγ1(ω(NL)1) and Tγ1α(NR); that is, MR and Wu(ω(NL)×(1ρ,1+ρ)) intersect transversally. We thus conclude that, there exists 0<ε0ε1 such that, if 0<εε0, then ML(ε,δ) intersects MR(ε) transversally.

    For uniqueness, note that the transversality of the intersection ML(ε,δ)MR(ε) implies dim(ML(ε,δ)MR(ε))=dimML(ε,δ)+dimMR(ε)9=1. Thus, there exists δ0>0 such that, if 0<δδ0, the intersection ML(ε,δ)MR(ε) consists of precisely one solution near the singular orbit Γ0ΛΓ1.

    In this section, our main interest is the competition between two positively charged ion species due to the nonlinear interplays between finite ion sizes, diffusion coefficients and boundary conditions. Mathematically, this is characterized by the quantity P=D1J11D2J21 (J11 and J21 are the leading terms corresponding to two cations that include finite size effects), which provides important information of competitions between two cations. For convenience, in the following, we use SI to denote the cation associated with c1 and J1, and SII to denote the one associated with c2 and J2.

    Recall from (3.7) that Jk1's (k=1,2) are the main terms that contain ion size effects. From (3.15), Lemmas 3.4 and 3.5, one has

    J11=[L3J10+J20(A1+Tc0Tm0(1)1)]1{(1+F2)J10L232(J10+J20)(A2+Tc0Tm0(1)1)+F1(Tc0Tm1Tc1Tm0)(Tm0)2lnA(1)+J10(Tm1Tc0Tc1Tm0)L32Tm0(J10+J20)2(A1+Tc0Tm0(1)1)+M1(1)M3(1)+(λ2λ1)F21(1+Tc02J30)(ATc0Tm0(1)1)},J21=[L3J10+J20(A1+Tc0Tm0(1)1)]1{(1+F2)J20L232(J10+J20)(A2+Tc0Tm0(1)1)+F1(Tc0Tm1Tc1Tm0)(Tm0)2lnA(1)+J20(Tm1Tc0Tc1Tm0)L32Tm0(J10+J20)2(A1+Tc0Tm0(1)1)+M2(1)+M4(1)(λ2λ1)F21(1+Tc02J30)(ATc0Tm0(1)1)}.

    For simplicity in our following discussion, we further assume that λ2>λ1. Define P(V)=P(L1,L2,R1,R2;V;λ1,λ2;D1,D2) as

    P=D1J11(L1,L2,R1,R2,V;λ1,λ2)D2J21(L1,L2,R1,R2,V;λ1,λ2).

    From Lemma 3.5, P(V) can be rewritten as

    P(V)=P0(V)R3eVL3, (4.1)

    with P0(V)=P0(L1,L2,R1,R2;V;λ1,λ2;D1,D2) being defined by

    P0=(D2J20D1J10)N1(D1J10+D2J20)(λ2λ1)F1Tc0L3Tm0(A(1)1)+(D1+D2)(J10+J20)N2, (4.2)

    where N1(V)=N1(L1,L2,R1,R2,V;λ1,λ2) and N2(V)=N2(L1,L2,R1,R2,V;λ1,λ2) are defined as follows:

    N1(V)=L3(ATc0Tm0(1)1)(Tm1Tm0(λ2λ1)F1Tc02J30L3(1+F2)2)+(1+F2)L232(A2+Tc0Tm0(1)1)(λ2λ1)L3F1(Tc0)22Tm0J30(A(1)1)(Tm1Tc0Tc1Tm0)L32Tm0(J10+J20)(A1+Tc0Tm0(1)1),N2(V)=F1Tc0Tm0(A1(1)1)(Tm1Tm0(λ2λ1)F1Tc02J30L3(1+F2)2)F1L3Tm0(1+F2)(Tm0+Tc02)(A(1)1)+F1(Tc0Tm1Tc1Tm0)(Tm0)2lnA(1)+(λ2λ1)F21(1+Tc02J30)(ATc0Tm0(1)1). (4.3)

    For the function P defined in (4.1), direct calculations lead to the following result, which is crucial for our further study on the qualitative properties of ionic flows.

    Lemma 4.1. Assume that V0. Under electroneutrality conditions, L3=R3, L3=R3eV and L3=R3eV are three removable singularities of P(V).

    Remark 4.2. From the definition of the function P(V), L3=R3, L3=R3eV and L3=R3eV are removable singularities of the first order approximations Jk1, k=1,2 in ν of the individual fluxes.

    We next consider the competition of two positive charged ion species as L3R3 and L3R3eV, respectively.

    Theorem 4.3. Suppose L1>R1. For small ε>0 and ν>0, there exists a critical potential V>ln2 such that, for V>V, one has P(V)>0 as L3R3, that is, the ion channel prefers cation SI over cation SII.

    Proof. Direct calculation gives

    limL3R3P(V)=(λ2λ1)(L1R1)VeV1[(D1+D2)(L1R1)1eV+D1R1+D2R2(D1L1+D2L2)eV(1eV)2V+(D2R3+(D1+D2)(R1L1eV)1eV)(eV12V1eV+(1eV)2V)].

    Under our assumption, one has (λ2λ1)(L1R1)VeV1<0, for V>0. Upon introducing K1(V) and K2(V) defined as

    K1(V)=eV12+(1eV)2V,K2(V)=eV12+(1eV)2VV1eV.

    Then, as L3R3, P(V)>0 is equivalent to

    D2V(R2R1(L2L1)eV)(1eV)2+D2R3K2(V)+(D1+D2)(L1R1)1eV+(D1+D2)(R1L1eV)1eVK1(V)<0.

    Rewrite K1(V) as K1(V)=eV12V˜K1(V) with ˜K1(V)=V+2eV2. A direct computation gives ˜K1(0)=0 and ˜K1(V)>0 for V>ln2; and ˜K1(V)<0 for V<ln2. Together with eV12V<0, one has K1(V)<0 for V>ln2 and K1(V)>0 for V<ln2. Note that K2(V)=K1(V)+VeV1 with VeV1>0. Therefore, there must exist a V>ln2 such that K2(V)<0 for V>V. With L1>R1 and V>0, we have

    (ⅰ) R1L1eV=R1L1+L1(1eV)<0, which implies (D1+D2)(R1L1eV)1eV>0.

    (ⅱ) R2R1<L2L1<(L2L1)eV gives D2V(1eV)2(R2R1(L2L1)eV)<0.

    (ⅲ) (D1+D2)(L1R1)1eV<0.

    Therefore, with L1>R1, one has, for V>V,

    D2V(1eV)2(R2R1(L2L1)eV)+D2R3K2(V)+(D1+D2)(L1R1)1eV+(D1+D2)(R1L1eV)1eVK1(V)<0.

    This completes the proof.

    Theorem 4.4. Suppose D1D2>L2L1 and V>0 (resp. D1D2<L2L1 and V<0). For small ε>0 and ν>0, one has P(V)>0 as L3R3eV, that is, the ion channel prefers cation SI over cation SII as L3R3eV.

    Proof. It follows from P(V)(D1L1D2L2)(eV1)2(λ1L1+λ2L2+R3eV)V as L3R3eV.

    For fixed boundary concentration R3 and boundary electric potential V>0, the removable singularities L3=R3, L3=R3eV and L3=R3eV split the L3-region into four subregions: (0,R3eV), (R3eV,R3), (R3,R3eV)  and  (R3eV,+), over which the sign of P can be analyzed. Our study will mainly focus on the sign of P over the first subregion (0,R3eV).

    To get started, we establish some results which are crucial to study the sign of P(V), and will be frequently used in our following analysis. Note that, under electroneutrality conditions (3.4), from (3.18), one has

    F1=J20L1J10L2J10+J20=L1R2L2R1R3L3eV, (4.4)

    and from (3.15), (5.3) and (5.4), together with A(1)=R3/L3, we have

    Tm1=((λ2λ1)F1Tc0J30(R3eVL3)+(1+F2)(R23L23)),Tc1=(λ2λ1)F1Tc0Tm02J30(lnR3lnL3)q1(V;L3,R3)+(1+F2)(L3+R3)Tc02q2(L3,R3), (4.5)

    where

    q1(V;L3,R3)=R3R3L3(lnR3lnL3)(R3L3)(lnR3lnL3)(eV1)V+eV1V,q2(L3,R3)=11L3R3L3lnR3lnL3.

    It then follows that

    Lemma 4.5. Under electroneutrality conditions, one has

    (ⅰ) For L3<R3eV, one has F1<0 (resp. F1>0) if L1L2<R1R2 (resp. L1L2>R1R2).

    (ⅱ) For L3>R3eV, one has F1<0 (resp. F1>0) if L1L2>R1R2 (resp. L1L2<R1R2).

    Lemma 4.6. Fixing V>0. Then,

    (Ⅰ) For R3<1, one has

    (ⅰ) q1(V;L3,R3)>0 if L3<R3;

    (ⅱ) There exists an L3>R3 such that q1(V;L3,R3)<0 if R3<L3<L3 and q1(V;L3,R3)>0 if L3>L3.

    (Ⅱ) For R3>1, one has q1(V;L3,R3)>0 if L3<R3; q1(V;L3,R3)<0 if L3>R3.

    Lemma 4.7. q2(L3,R3)>0 if L3>R3; and q2(L3,R3)<0 if L3<R3.

    Based on Lemma 4.6, we further assume R3>1 in the following argument (the case with R3<1 can be discussed similarly). From the electroneaurality conditions (3.4), one has the following possibilities for L3<R3eV:(a) L1<R1eV and L2<R2eV; (b) L1<R1eV and L2>R2eV; and (c) L1>R1eV and L2<R2eV. In this section, we will only consider the case (a), which, with V>0, is equivalent to

    0<V<min{lnR1L1, lnR2L2}, (4.6)

    and study the sign of P, which provides rich dynamical behavior of ionic flows. The other two cases can be analyzed in a similar way, and we leave them to interested readers. From Lemmas 4.5–4.7, we obtain

    Lemma 4.8. Under assumption (4.6), one has J10<0,  J20<0,  J30<0,  Tm0<0 and Tc0>0. Furthermore, if L1L2>R1R2, then Tm1<0 and Tc1<0.

    Our main result now follows.

    Theorem 4.9. Assume (4.6), L1L2>R1R2 and R1L1R2L2<D2D1. Then, for small ε>0 and ν>0, there exists a critical potential V>0 such that for 0<V<V, one has P(V)>0, that is, the ion channel prefers cation SI over cation SII.

    Proof. Note that, under our assumption, PP0>0. Instead of showing that P>0, we will prove that P0>0. Recall from (4.2),

    P0=(D2J20D1J10)N1(D1J10+D2J20)(λ2λ1)F1Tc0L3Tm0(A(1)1)+(D1+D2)(J10+J20)N2,

    where N1 and N2 as defined in (4.3). From Lemma 4.8, one can easily check that for the second term in P0, one has (D1J10+D2J20)(λ2λ1)F1Tc0L3Tm0(A(1)1)>0. Next we show that, under our assumption, (D2J20D1J10)N1+(D1+D2)(J10+J20)N2>0 for 0<V<V, where V is identified at the end of our proof. For convenience, we rewrite N1 and N2 as

    N1=L3(ATc0Tm0(1)1)(Tm1Tm0(λ2λ1)F1Tc02J30L3(1+F2)2)+N10,N2=F1Tc0Tm0(A1(1)1)(Tm1Tm0(λ2λ1)F1Tc02J30L3(1+F2)2)+N20,

    where

    N10=(1+F2)L232(A2+Tc0Tm0(1)1)(Tm1Tc0Tc1Tm0)L32Tm0(J10+J20)(A1+Tc0Tm0(1)1)(λ2λ1)L3F1(Tc0)22Tm0J30(A(1)1),N20=(λ2λ1)F21(1+Tc02J30)(ATc0Tm0(1)1)+F1(Tc0Tm1Tc1Tm0)(Tm0)2lnA(1)F1L3Tm0(1+F2)(Tm0+Tc02)(A(1)1).

    Notice that, with L1<R1eV and L2<R2eV,

    A2+Tc0Tm0(1)1=R3L3R3L3eV1>0,A1+Tc0Tm0(1)1=R3L3eV1>0,ATc0Tm0(1)1=eV1>0,A(1)1=R3L31>0,Tm0+Tc02=Tm0(1+V2(lnL3lnR3))<0,1+Tc02J30=lnL3lnR3lnL3lnR3V>0.

    Together with Lemma 4.8, one has N10>0  and  N20>0.

    Finally we prove that

    [(D2J20D1J10)L3(ATc0Tm0(1)1)+(D1+D2)(J10+J20)F1Tc0Tm0(A1(1)1)]×(Tm1Tm0(λ2λ1)F1Tc02J30L3(1+F2)2)>0.

    From (3.15), (3.18) and (5.3),

    Tm1Tm0(λ2λ1)F1Tc02J30L3(1+F2)2=(λ2λ1)(eV1)R3Tc0F12J30(R3L3)+(1+F2)R32>01+F2>(λ2λ1)Tc0(1eV)F1J30(R3L3),

    which is true under our assumption. It suffices to show that

    (D2J20D1J10)L3(ATc0Tm0(1)1)+(D1+D2)(J10+J20)F1Tc0Tm0(A1(1)1)>0,

    which is, from Lemma 3.4 and Lemma 4.8, equivalent to ω(V)<0, where ω(V):=ω(V;L1,L2,L3,R1,R2,R3;D1,D2) is defined as

    ω(V)=(D2(R2L2eV)D1(R1L1eV))L3(eV1)+(D1+D2)(L1R2L2R1)(L3R3)V(lnL3lnR3)R3.

    With R1L1R2L2<D2D1, a careful calculation leads to ω(0)=0; ω(0)<0 and ω(V)>0 for V>0. It follows that the curve ω(V) is concave up for V>0 with its minimum being negative. Therefore, there exists a unique Vg such that ω(V)<0 for 0<V<Vg. Define V=min{Vg,lnR1L1,lnR2L2}, one has P(V)>0 for 0<V<V.

    Remark 4.10. In Theorem 4.9 (see also Theorems 4.3 and 4.4 in Subsection 4.1), under electroneutrality boundary concentration conditions, we study the competition between two cations in terms of the function P defined in (4.1), which provides useful insights for the study of ionic flows through membrane channels, in particular, the selectivity of cations due to the electrodiffusion property of ions with given protein structures of ion channels (another important factor in studying ion selectivity). Distinct effects of the nonlinear interplays between the physical parameters, such as relative ion size (λ1,λ2), diffusion coefficients (D1,D2), boundary concentrations (L1,L2,L3), (R1,R2,R3) and boundary potential V are characterized in a great detail. For example, under the assumption of (4.6), L1L2<R1R2, and R1L1R2L2<D2D1, one is able to identify some critical potential V, such that for 0<V<V, the ion channel prefers cation SI over cation SII. Together with the scaling law stated in Proposition 3.6 and further illustrated in Remark 3.7, the preference of ion channels for distinct cations could be either reduced or enhanced by choosing some suitable positive parameter α in Proposition 3.6.

    To end this section, we perform numerical simulations for the system (3.1) directly with small ε and ν to detect two critical potentials identified in Theorem 4.3 and Theorem 4.9, respectively (Figures 2 and 3).

    Figure 2.  Numerical simulations to detect the critical potential V identified in Theorem 4.3 with ε=0.008 and ν=0.001 as L3R3 (left figure) and L3R+3 (right one), respectively.
    Figure 3.  Numerical simulations to detect the critical potential V identified in Theorem 4.9 with ε=0.008 and ν=0.001. Clearly, V<min{lnR1L1=1.6964, lnR2L2=1.7918}.

    We studied a quasi-one-dimensional PNP model with three ion species, two positively charged with the same valence and one negatively charged, and Bikerman's local hard-sphere potential accounted for ion size effects on ionic flows. Under the framework of a geometric singular perturbation theory, the existence and uniqueness of the boundary value problem was established. Furthermore, treating the ion sizes as small parameters, both zeroth order and first order approximations (in ν) to individual fluxes are derived. Of particular interest is to examine the competition between two positively charged ion species for open ion channels with given protein structures, which provides useful insights into the study of selectivity phenomena. Based on our rigorous analysis, we are able to characterize the distinct effects of the nonlinear interplay between physical parameters, such as relative ion size (λ1,λ2), diffusion coefficients (D1,D2), boundary concentrations (L1,L2,L3), (R1,R2,R3) and boundary potential V, which provides an efficient way to control ionic flows. These are the novelty and main contribution of our work. The results, although established for simple biology settings have demonstrated extremely rich dynamics of ionic flows and sensitive dependence on all those physical parameters. More complex phenomena for more realistic ion channel models are expected. We would like to point out that the topic studied in this work is closely related to artificial nanopores as well, where modeling mainly focuses on computing device functions: ionic currents as responses to the voltage as a stimulus (see [36,38,78,89,101] for more details). We believe this work will be useful for future numerical studies, stimulate further analytical studies of ionic flows concerning the selectivity of cations, and even provide useful insights into biological experiments.

    Appendix A. We provide detailed proof for Lemma 3.5 in section 3.2.

    Proof of Lemma 3.5 From c11 and c21 equations in (3.17), one has

    ˙c11+˙c21=(λ1+1)Tm0+2(λ2λ1)J202c10+(λ2+1)Tm0+2(λ1λ2)J102c20Tm12.

    Together with Proposition 3.3, one has

    c11(x)+c21(x)=Tm12x+(λ1+1)Tm0+2(λ2λ1)J202x0c10(s)ds+(λ2+1)Tm0+2(λ1λ2)J102x0c20(s)ds.

    From (3.18), together with x=2L3Tm0(1A(x)), we have

    x0c10(s)ds=F1L3J30(A1Tc0Tm0(x)1)(A2(x)1)L23Tm0(J10+J20)J10,x0c20(s)ds=F1L3J30(A1Tc0Tm0(x)1)(A2(x)1)L23Tm0(J10+J20)J20. (5.1)

    Therefore,

    c11(x)+c21(x)=Tm12x(λ2λ1)F1Tc0L32J30(A1Tc0Tm0(x)1)L232(1+F2)(A2(x)1). (5.2)

    Evaluating Eq (5.2) at x=1, we have

    0=Tm12(λ2λ1)F1Tc0L32J30(A1Tc0Tm0(1)1)L232(1+F2)(A2(1)1).

    It follows that

    Tm1=L3[(λ2λ1)F1Tc0J30(A1Tc0Tm0(1)1)+L3(1+F2)(A2(1)1)]. (5.3)

    Solving the ϕ1 equations in (3.17), we have

    ϕ1(x)=Tc02x0c11(s)+c21(s)(c10(s)+c20(s))2dsTc12x01c10(s)+c20(s)ds.

    Note that, from (3.18) and (5.2), together with 1A(x)=Tm02L3x, one has

    x0dsc10(s)+c20(s)=2Tm0lnA(x),x0(c11(s)+c21(s))ds(c10(s)+c20(s))2=L3(1+F2)Tm0(A(x)1)(λ2λ1)F1J30(ATc0Tm0(x)1)2Tm1(Tm0)2lnA(x)+2Tm0(L3(1+F2)2+(λ2λ1)Tc0F12J30Tm1Tm0)(A1(x)1).

    Therefore,

    ϕ1(x)=Tc0Tm0(L3(1+F2)2+(λ2λ1)F1Tc02J30Tm1Tm0)(A1(x)1)+L3(1+F2)I0(A(x)1)2Tm0(λ2λ1)F1Tc0(ATc0Tm0(x)1)2J30+Tc1Tm0Tm1Tc0(Tm0)2lnA(x).

    Evaluating ϕ1(x) at x=1 gives

    Tc1=Tm0lnA(1)[Tc0Tm0(Tm1Tm0L3(1+F2)2(λ2λ1)F1Tc02J30)(A1(1)1)L3(1+F2)Tc02Tm0(A(1)1)+(λ2λ1)F1Tc02J30(ATc0Tm0(1)1)]+Tm1Tc0Tm0. (5.4)

    By the variation of constant formula, one has, with cL11=cL21=0,

    c11(x)=ATc0Tm0(x)[Tc02x0c10(s)(c11(s)+c21(s))ATc0Tm0(s)(c10(s)+c20(s))2dsJ11x0ATc0Tm0(s)ds+Tc12x0c10(s)ATc0Tm0(s)dsc10(s)+c20(s)+(λ1+1)Tm0+2(λ2λ1)J202x0c210(s)ATc0Tm0(s)dsc10(s)+c20(s)+(λ2+1)Tm0+2(λ1λ2)J102x0c10(s)c20(s)ATc0Tm0(s)dsc10(s)+c20(s)],c21(x)=ATc0Tm0(x)[Tc02x0c20(s)(c11(s)+c21(s))ATc0Tm0(s)(c10(s)+c20(s))2dsJ21x0ATc0Tm0(s)ds+Tc12x0c20(s)ATc0Tm0(s)dsc10(s)+c20(s)+(λ2+1)Tm0+2(λ1λ2)J102x0c220(s)ATc0Tm0(s)dsc10(s)+c20(s)+(λ1+1)Tm0+2(λ2λ1)J202x0c10(s)c20(s)ATc0Tm0(s)dsc10(s)+c20(s)],

    where, from (3.18), (5.1) and (5.2),

    x0c10(s)(c11(s)+c21(s))(c10(s)+c20(s))2ds=2Tm0{J10Tm1L32(J10+J20)2(A1+Tc0Tm0(x)1)+Tm0Tc0M1(x)+(λ2λ1)Tm0F212J30(ATc0Tm0(x)1)(1+F2)J10Tm0L232(J10+J20)(Tc0+2Tm0)(A2+Tc0Tm0(x)1)+F1Tm1Tm0lnA(x)L3F12((λ2λ1)J10Tc0(J10+J20)J30+1+F2)(A(x)1)},x0c10(s)ATc0Tm0(s)c10(s)+c20(s)ds=2F1Tm0lnA(x)J10L3(J10+J20)2(A1+Tc0Tm0(x)1),x0c210(s)ATc0Tm0(s)dsc10(s)+c20(s)=2Tm0[(J10L3J10+J20)2Tm0Tc0+2Tm0(A2+Tc0Tm0(x)1)+2F1J10L3J10+J20(A(x)1)F21Tm0I0(ATc0Tm0(x)1)],x0c20(s)(c11(s)+c21(s))(c10(s)+c20(s))2ds=2Tm0{J20Tm1L32(J10+J20)2(A1+Tc0Tm0(x)1)+Tm0Tc0M2(x)(λ2λ1)Tm0F212J30(ATc0Tm0(x)1)(1+F2)J20Tm0L232(J10+J20)(Tc0+2Tm0)(A2+Tc0Tm0(x)1)F1Tm1Tm0lnA(x)+L3F12(1+F2(λ2λ1)J20Tc0(J10+J20)J30)(A(x)1)},x0c20(s)ATc0Tm0(s)c10(s)+c20(s)ds=2F1Tm0lnA(x)J20L3(J10+J20)2(A1+Tc0Tm0(x)1),x0c220(s)ATc0Tm0(s)dsc10(s)+c20(s)=2Tm0[(J20L3J10+J20)2Tm0Tc0+2Tm0(A2+Tc0Tm0(x)1)2F1J20L3J10+J20(A(x)1)F21Tm0Tc0(ATc0Tm0(x)1)],x0c10(s)c20(s)ATc0Tm0(s)dsc10(s)+c20(s)=2Tm0[J10J20Tm0L23(J10+J20)2(Tc0+2Tm0)(A2+Tc0Tm0(x)1)+F21Tm0I0(ATc0Tm0(x)1)+F1L3(J20J10)J10+J20(A(x)1)].

    Evaluating ϕ1(x), c11(x) and c21(x) at x=1, the expressions for J11,J21 and J31 follows directly.

    The authors are grateful to the anonymous referees whose suggestions have in our opinion, significantly improved the paper. J. Chen and M. Zhang were partially supported by MPS Simons Foundation (No. 628308).

    The author declares no conflicts of interest in this paper.



    [1] N. Abaid, R. S. Eisenberg, W. Liu, Asymptotic expansions of I-V relations via a Poisson-NernstPlanck system, SIAM J. Appl. Dyn. Syst., 7 (2008), 1507-1526. doi: 10.1137/070691322
    [2] R. Aitbayev, P. W. Bates, H. Lu, L. Zhang, M. Zhang, Mathematical studies of Poisson-NernstPlanck systems: dynamics of ionic flows without electroneutrality conditions, J. Comput. Appl. Math., 362 (2019), 510-527. doi: 10.1016/j.cam.2018.10.037
    [3] V. Barcilon, Ion flow through narrow membrane channels: Part I, SIAM J. Appl. Math., 52 (1992), 1391-1404. doi: 10.1137/0152080
    [4] V. Barcilon, D. P. Chen, R. S. Eisenberg, Ion flow through narrow membrane channels: Part II, SIAM J. Appl. Math., 52 (1992), 1405-1425. doi: 10.1137/0152081
    [5] V. Barcilon, D. P. Chen, R. S. Eisenberg, J. W. Jerome, Qualitative properties of steady-state Poisson-Nernst-Planck systems: Perturbation and simulation study, SIAM J. Appl. Math., 57 (1997), 631-648. doi: 10.1137/S0036139995312149
    [6] P. W. Bates, Y. Jia, G. Lin, H. Lu, M. Zhang, Individual flux study via steady-state PoissonNernst-Planck systems: Effects from boundary conditions, SIAM J. Appl. Dyn. Syst., 16 (2017),
    [7] P. W. Bates, W. Liu, H. Lu, M. Zhang, Ion size and valence effects on ionic flows via PoissonNernst-Planck models, Commu. Math. Sci., 15 (2017), 881-901. doi: 10.4310/CMS.2017.v15.n4.a1
    [8] J. J. Bikerman, Structure and capacity of the electrical double layer, Philos. Mag. J. Sci., 33 (1942), 384-397. doi: 10.1080/14786444208520813
    [9] D. Boda, D. Busath, B. Eisenberg, D. Henderson, W. Nonner, Monte Carlo simulations of ion selectivity in a biological Na+ channel: Charge-space competition, Phys. Chem. Chem. Phys., 4 (2002), 5154-5160. doi: 10.1039/B203686J
    [10] D. Boda, D. Gillespie, W. Nonner, D. Henderson, B. Eisenberg, Computing induced charges in inhomogeneous dielectric media: Application in a Monte Carlo simulation of complex ionic systems, Phys. Rev. E, 69 (2004), 046702.
    [11] D. Boda, W. Nonner, M. Valiskó, D. Henderson, B. Eisenberg, D. Gillespie, Steric selectivity in Na+ channels arising from protein polarization and mobile side chains, Biophys J., 93 (2007), 1960-1980. doi: 10.1529/biophysj.107.105478
    [12] D. Boda, W. Nonner, D. Henderson, B. Eisenberg, D. Gillespie, Volume exclusion in calcium selective channels, Biophys J., 94 (2008), 3486-3496. doi: 10.1529/biophysj.107.122796
    [13] M. Burger, R. S. Eisenberg, H. W. Engl, Inverse problems related to ion channel selectivity, SIAM J. Appl. Math., 67 (2007), 960-989. doi: 10.1137/060664689
    [14] A. E. Cardenas, R. D. Coalson, M. G. Kurnikova, Three-Dimensional Poisson-Nernst-Planck Theory Studies: Influence of Membrane Electrostatics on Gramicidin A Channel Conductance, Biophys. J., 79 (2000), 80-93. doi: 10.1016/S0006-3495(00)76275-8
    [15] J. H. Chaudhry, S. D. Bond, L. N. Olson, Finite Element Approximation to a Finite-Size Modified Poisson-Boltzmann Equation, J. Sci. Comput., 47 (2011), 347-364. doi: 10.1007/s10915-010-9441-7
    [16] D. P. Chen, R. S. Eisenberg, Charges, currents and potentials in ionic channels of one conformation, Biophys. J., 64 (1993), 1405-1421.
    [17] S. Chung, S. Kuyucak, Predicting channel function from channel structure using Brownian dynamics simulations, Clin. Exp. Pharmacol Physiol., 28 (2001), 89-94. doi: 10.1046/j.1440-1681.2001.03408.x
    [18] J. R. Clay, Potassium current in the squid giant axon, Int. Rev. Neurobiol., 27 (1985), 363-384. doi: 10.1016/S0074-7742(08)60562-0
    [19] R. Coalson, M. Kurnikova, Poisson-Nernst-Planck theory approach to the calculation of current through biological ion channels, IEEE Trans. NanoBioscience, 4 (2005), 81-93. doi: 10.1109/TNB.2004.842495
    [20] B. Corry, T. W. Allen, S. Kuyucak, S. H. Chung, Mechanisms of permeation and selectivity in calcium channels, Biophys J., 80 (2001), 195-214. doi: 10.1016/S0006-3495(01)76007-9
    [21] B. Corry, T. W. Allen, S. Kuyucak, S. H. Chung, A model of calcium channels, Biochim. Biophys. Acta Biomembr., 1509 (2000), 1-6.
    [22] B. Corry, S. H. Chung, Mechanisms of valence selectivity in biological ion channels, Cell. Mol. Life Sci., 63 (2006), 301-315. doi: 10.1007/s00018-005-5405-8
    [23] J. M. Diamond, E. M. Wright, Biological membranes: the physical basis of ion and nonelectrolyte selectivity, Annu. Rev. Physiol., 31 (1969), 581-646. doi: 10.1146/annurev.ph.31.030169.003053
    [24] D. A. Doyle, J. M. Cabral, R. A. Pfuetzner, A. Kuo, J. M. Gulbis, S. L. Cohen, et al., The Structure of the Potassium Channel: Molecular Basis of K+ Conduction and Selectivity, Science, 280 (1998), 69-77.
    [25] R. Dutzler, E. B. Campbell, M. Cadene, B. T. Chait, R. Mackinnon, X-ray structure of a ClC chloride channel at 3.0 A reveals the molecular basis of anion selectivity, Nature, 415 (2002), 287-294. doi: 10.1038/415287a
    [26] R. Dutzler, E.B. Campbell, R. MacKinnon, Gating the selectivity filter in ClC chloride channels, Science, 300 (2003), 108-112. doi: 10.1126/science.1082708
    [27] B. Eisenberg, Ion Channels as Devices, J. Comput. Electron., 2 (2003), 245-249.
    [28] B. Eisenberg, Proteins, Channels, and Crowded Ions, Biophys. Chem., 100 (2003), 507-517.
    [29] R. S. Eisenberg, Channels as enzymes, J. Memb. Biol., 115 (1990), 1-12.
    [30] R. S. Eisenberg, R. Elber, Atomic biology, electrostatics and Ionic Channels, in New Developments and Theoretical Studies of Proteins, World Scientific, (1996), 269-357.
    [31] R. S. Eisenberg, From Structure to Function in Open Ionic Channels, J. Memb. Biol., 171 (1999), 1-24. doi: 10.1007/s002329900554
    [32] B. Eisenberg, W. Liu, Poisson-Nernst-Planck systems for ion channels with permanent charges, SIAM J. Math. Anal., 38 (2007), 1932-1966. doi: 10.1137/060657480
    [33] G. Eisenman, R. Horn, Ionic selectivity revisited: The role of kinetic and equilibrium processes in ion permeation through channels, J. Memb. Biol., 76 (1983), 197-225. doi: 10.1007/BF01870364
    [34] A. Ern, R. Joubaud, T. Leliévre, Mathematical study of non-ideal electrostatic correlations in equilibrium electrolytes, Nonlinearity, 25 (2012), 1635-1652. doi: 10.1088/0951-7715/25/6/1635
    [35] N. Fenichel, Geometric singular perturbation theory for ordinary differential equations, J. Differ. Equations, 31 (1979), 53-98. doi: 10.1016/0022-0396(79)90152-9
    [36] D. Fertig, B. Matejczyk, M. Valiskó, D. Gillespie, D. Boda, Scaling Behavior of Bipolar Nanopore Rectification with Multivalent Ions, J. Phys. Chem. C., 123 (2019), 28985-28996. doi: 10.1021/acs.jpcc.9b07574
    [37] J. Fischer, U. Heinbuch, Relationship between free energy density functional, Born-Green-Yvon, and potential distribution approaches for inhomogeneous fluids, J. Chem. Phys., 88 (1988), 1909-1913.
    [38] T. Gamble, K. Decker, T. S Plett, M. Pevarnik, J.F. Pietschmann, I. V. Vlassiouk, et al., Rectification of ion current in nanopores depends on the type of monovalent cations-experiments and modeling, J. Phys. Chem. C, 118 (2014), 9809-9819.
    [39] D. Gillespie, A singular perturbation analysis of the Poisson-Nernst-Planck system: Applications to Ionic Channels, Ph. D Dissertation, Rush University, Chicago, 1999.
    [40] D. Gillespie, R. S. Eisenberg, Physical descriptions of experimental selectivity measurements in ion channels, European Biophys. J., 31 (2002), 454-466. doi: 10.1007/s00249-002-0239-x
    [41] D. Gillespie, W. Nonner, R. S. Eisenberg, Coupling Poisson-Nernst-Planck and density functional theory to calculate ion flux, J. Phys. Condens. Matter, 14 (2002), 12129-12145. doi: 10.1088/0953-8984/14/46/317
    [42] D. Gillespie, W. Nonner, R. S. Eisenberg, Density functional theory of charged, hard-sphere fluids, Phys. Rev. E, 68 (2003), 0313503.
    [43] D. Gillespie, W. Nonner, R. S. Eisenberg, Crowded Charge in Biological Ion Channels, Nanotech, 3 (2003), 435-438.
    [44] E. Gongadze, U. van Rienen, V. Kralj-lglič, A. lglič, Spatial variation of permittivity of an electrolyte solution in contact with a charged metal surface: A mini review, Comput. Method Biomech. Biomed. Eng., 16 (2013), 463-480.
    [45] B. Hille, Ionic Channels of Excitable Membranes, Sinauer Associates Inc, (2001).
    [46] B. Hille, Ionic Selectivity, saturation, and block in sodium channels. A four barrier model, J. Gen. Physiol., 66 (1975), 535-560.
    [47] U. Hollerbach, D. P. Chen, R. S. Eisenberg, Two- and Three-Dimensional Poisson-Nernst-Planck Simulations of Current Flow through Gramicidin-A, J. Comp. Science, 16 (2002), 373-409.
    [48] U. Hollerbach, D. Chen, W. Nonner, B. Eisenberg, Three-dimensional Poisson-Nernst-Planck Theory of Open Channels, Biophys. J., 76 (1999), A205.
    [49] A. L. Hodgkin, The Conduction of the Nervous Impulse, Liverpool University Press, (1971), 1-108.
    [50] A. L. Hodgkin, Chance and Design, Cambridge University Press, (1992).
    [51] A. L. Hodgkin, A. F. Huxley, Propagation of electrical signals along giant nerve fibres, Proc. R. Soc. London B, 140 (1952), 177-183. doi: 10.1098/rspb.1952.0054
    [52] A. L. Hodgkin, A. F. Huxley, Currents carried by sodium and potassium ions through the membrane of the giant axon of Loligo, J. Physiol., 116 (1952), 449-472. doi: 10.1113/jphysiol.1952.sp004717
    [53] A. L. Hodgkin, A. F. Huxley, The components of membrane conductance in the giant axon of Loligo, J. Physiol., 116 (1952), 473-496. doi: 10.1113/jphysiol.1952.sp004718
    [54] A. L. Hodgkin, A. F. Huxley, The dual effect of membrane potential on sodium conductance in the giant axon of Loligo, J. Physiol., 116 (1952), 497-506. doi: 10.1113/jphysiol.1952.sp004719
    [55] A. L. Hodgkin, A. F. Huxley, A quantitive description of membrane current and its application to conduction and excitation in nerve, J. Physiol., 117 (1952), 500-544. doi: 10.1113/jphysiol.1952.sp004764
    [56] Y. Hyon, B. Eisenberg, C. Liu, A mathematical model for the hard sphere repulsion in ionic solutions, Commun. Math. Sci., 9 (2010), 459-475.
    [57] Y. Hyon, J. Fonseca, B. Eisenberg, C. Liu, Energy variational approach to study charge inversion (layering) near charged walls, Discrete Contin. Dyn. Syst. B, 17 (2012), 2725-2743. doi: 10.3934/dcdsb.2012.17.2725
    [58] Y. Hyon, C. Liu, B. Eisenberg, PNP equations with steric effects: a model of ion flow through channels, J. Phys. Chem. B, 116 (2012), 11422-11441. doi: 10.1021/jp305273n
    [59] W. Im, D. Beglov, B. Roux, Continuum solvation model: Electrostatic forces from numerical solutions to the Poisson-Bolztmann equation, Comp. Phys. Comm., 111 (1998), 59-75. doi: 10.1016/S0010-4655(98)00016-2
    [60] W. Im, B. Roux, Ion permeation and selectivity of OmpF porin: A theoretical study based on molecular dynamics, Brownian dynamics, and continuum electrodiffusion theory, J. Mol. Biol., 322 (2002), 851-869.
    [61] S. Ji, W. Liu, Poisson-Nernst-Planck Systems for Ion Flow with Density Functional Theory for Hard-Sphere Potential: I-V relations and Critical Potentials. Part I: Analysis, J. Dyn. Diff. Equat., 24 (2012), 955-983. doi: 10.1007/s10884-012-9277-y
    [62] S. Ji, W. Liu, M. Zhang, Effects of (small) permanent charges and channel geometry on ionic flows via classical Poisson-Nernst-Planck models, SIAM J. on Appl. Math., 75 (2015), 114-135. doi: 10.1137/140992527
    [63] Y. Jia, W. Liu, M. Zhang, Qualitative properties of ionic flows via Poisson-Nernst-Planck systems with Bikerman's local hard-sphere potential: Ion size effects, Discrete Contin. Dyn. Syst. B, 21 (2016), 1775-1802. doi: 10.3934/dcdsb.2016022
    [64] C. Jones, Geometric singular perturbation theory, in Dynamical systems, Springer, (1995), 44-118.
    [65] C. Jones, T. Kaper, N. Kopell, Tracking invariant manifolds up to exponentially small errors, SIAM J. Math. Anal., 27 (1996), 558-577. doi: 10.1137/S003614109325966X
    [66] C. Jones, N. Kopell, Tracking invariant manifolds with differential forms in singularly perturbed systems, J. Differ. Equations, 108 (1994), 64-88. doi: 10.1006/jdeq.1994.1025
    [67] A. S. Khair, T. M. Squires, Ion steric effects on electrophoresis of a colloidal particle, J. Fluid Mech., 640 (2009), 343-356. doi: 10.1017/S0022112009991728
    [68] M. S. Kilic, M. Z. Bazant, A. Ajdari, Steric effects in the dynamics of electrolytes at large applied voltages. II. Modified Poisson-Nernst-Planck equations, Phys. Rev. E, 75 (2007), 021503.
    [69] C. C. Lee, H. Lee, Y. Hyon, T. C. Lin, C. Liu, New Poisson-Boltzmann type equations: Onedimensional solutions, Nonlinearity, 24 (2011), 431-458. doi: 10.1088/0951-7715/24/2/004
    [70] B. Li, Continuum electrostatics for ionic solutions with non-uniform ionic sizes, Nonlinearity, 22 (2009), 811-833. doi: 10.1088/0951-7715/22/4/007
    [71] G. Lin, W. Liu, Y. Yi, M. Zhang, Poisson-Nernst-Planck systems for ion flow with density functional theory for local hard-sphere potential, SIAM J. Appl. Dyn. Syst., 12 (2013), 1613-1648. doi: 10.1137/120904056
    [72] W. Liu, Geometric singular perturbation approach to steady-state Poisson-Nernst-Planck systems, SIAM J. Appl. Math., 65 (2005), 754-766. doi: 10.1137/S0036139903420931
    [73] W. Liu, One-dimensional steady-state Poisson-Nernst-Planck systems for ion channels with multiple ion species, J. Differ. Equations, 246 (2009), 428-451. doi: 10.1016/j.jde.2008.09.010
    [74] W. Liu, H. Xu, A complete analysis of a classical Poisson-Nernst-Planck model for ionic flow, J. Differ. Equations, 258 (2015), 1192-1228. doi: 10.1016/j.jde.2014.10.015
    [75] W. Liu, B. Wang, Poisson-Nernst-Planck systems for narrow tubular-like membrane channels, J. Dyn. Differ. Equations, 22 (2010), 413-437. doi: 10.1007/s10884-010-9186-x
    [76] W. Liu, X. Tu, M. Zhang, Poisson-Nernst-Planck Systems for Ion Flow with Density Functional Theory for Hard-Sphere Potential: I-V relations and Critical Potentials. Part II: Numerics, J. Dyn. Differ. Equations, 24 (2012), 985-1004. doi: 10.1007/s10884-012-9278-x
    [77] H. Lu, J. Li, J. Shackelford, J. Vorenberg, M. Zhang, Ion size effects on individual fluxes via Poisson-Nernst-Planck systems with Bikerman's local hard-sphere potential: Analysis without electroneutrality boundary conditions, Discrete Contin. Dyn. Syst. B, 23 (2018), 1623-1643.
    [78] B. Matejczyk, M. Valiskó, M. T. Wolfram, J. F. Pietschmann, D. Boda, Multiscale modeling of a rectifying bipolar nanopore: Comparing Poisson-Nernst-Planck to Monte Carlo, J. Chem. Phys., 146 (2017), 124125.
    [79] J. H. Morais-Cabral, Y. Zhou, R. MacKinnon, Energetic optimization of ion conduction rate by the K+ selectivity filter, Nature, 414 (2001), 37-42. doi: 10.1038/35102000
    [80] L. J. Mullins, Ion selectivity of carriers and channels, Biophys J., 15 (1975), 921-931. doi: 10.1016/S0006-3495(75)85867-X
    [81] B. Nadler, Z. Schuss, A. Singer, B. Eisenberg, Diffusion through protein channels: from molecular description to continuum equations, Technical Proceedings of the Nanotechnology Conference and Trade Show, 2003. Available from: https: //ftp.rush.edu/.
    [82] W. Nonner, R. S. Eisenberg, Ion permeation and glutamate residues linked by Poisson-NernstPlanck theory in L-type Calcium channels, Biophys. J., 75 (1998), 1287-1305. doi: 10.1016/S0006-3495(98)74048-2
    [83] S. Y. Noskov, W. Im, B. Roux, Ion Permeation through the z1-Hemolysin Channel: Theoretical Studies Based on Brownian Dynamics and Poisson-Nernst-Planck Electrodiffusion Theory, Biophys. J., 87 (2004), 2299-2309. doi: 10.1529/biophysj.104.044008
    [84] S. Y. Noskov, B. Roux, Ion selectivity in potassium channels, Biophys. Chem., 124 (2006), 279-291.
    [85] J. K. Park, J. W. Jerome, Qualitative properties of steady-state Poisson-Nernst-Planck systems: Mathematical study, SIAM J. Appl. Math., 57 (1997), 609-630. doi: 10.1137/S0036139995279809
    [86] J. K. Percus, Equilibrium state of a classical fluid of hard rods in an external field, J. Stat. Phys., 15 (1976), 505-511. doi: 10.1007/BF01020803
    [87] J. K. Percus, Model grand potential for a nonuniform classical fluid, J. Chem. Phys., 75 (1981), 1316-1319. doi: 10.1063/1.442136
    [88] Y. Qiao, X. Liu, M. Chen, B. Lu, A Local Approximation of Fundamental Measure Theory Incorporated into Three Dimensional Poisson-Nernst-Planck Equations to Account for Hard Sphere Repulsion Among Ions, J. Stat. Phys., 163 (2016), 156-174. doi: 10.1007/s10955-016-1470-7
    [89] P. Ramirez, J. A. Manzanares, J. Cervera, V. Gomez, M. Ali, I. Pause, et al., Nanopore charge inversion and current-voltage curves in mixtures of asymmetric electrolytes, J. Membr. Sci., 563 (2018), 633-642.
    [90] Y. Rosenfeld, Free-Energy Model for the Inhomogeneous Hard-Sphere Fluid Mixture and Density-Functional Theory of Freezing, Phys. Rev. Lett., 63 (1989), 980-983. doi: 10.1103/PhysRevLett.63.980
    [91] Y. Rosenfeld, Free energy model for the inhomogeneous fluid mixtures: Yukawa-charged hard spheres, general interactions, and plasmas, J. Chem. Phys., 98 (1993), 8126-8148.
    [92] R. Roth, Fundamental measure theory for hard-sphere mixtures: A review, J. Phys. Condens. Matter, 22 (2010), 063102.
    [93] B. Roux, T. W. Allen, S. Berneche, W. Im, Theoretical and computational models of biological ion channels, Quat. Rev. Biophys., 37 (2004), 15-103. doi: 10.1017/S0033583504003968
    [94] D. Colquhoun, F. J. Sigworth, Single-Channel Recording, Plenum Press, (1995).
    [95] Z. Schuss, B. Nadler, R. S. Eisenberg, Derivation of Poisson and Nernst-Planck equations in a bath and channel from a molecular model, Phys. Rev. E, 64 (2001), 1-14.
    [96] F. Siddiqua, Z. Wang, S. Zhou, A Modified Poisson-Nernst-Planck Model with Excluded Volume Effect: Theory and Numerical Implementation, Commun. Math. Sci., 16 (2018), 251-271. doi: 10.4310/CMS.2018.v16.n1.a12
    [97] A. Singer, D. Gillespie, J. Norbury, R. S. Eisenberg, Singular perturbation analysis of the steadystate Poisson-Nernst-Planck system: applications to ion channels, Eur. J. Appl. Math., 19 (2008), 541-560. doi: 10.1017/S0956792508007596
    [98] A. Singer, J. Norbury, A Poisson-Nernst-Planck model for biological ion channels-an asymptotic analysis in a three-dimensional narrow funnel, SIAM J. Appl. Math., 70 (2009), 949-968. doi: 10.1137/070687037
    [99] L. Sun, W. Liu, Non-localness of Excess Potentials and Boundary Value Problems of PoissonNernst-Planck Systems for Ionic Flow: A Case Study, J. Dyn. Differ. Equations, 30 (2018), 779-797. doi: 10.1007/s10884-017-9578-2
    [100] S. K. Tin, N. Kopell, C. Jones, Invariant manifolds and singularly perturbed boundary value problems, SIAM J. Numer. Anal., 31 (1994), 1558-1576. doi: 10.1137/0731081
    [101] M. Valiskó, B. Matejczyk, Z. Ható, T. Kristóf, E. Mádai, D. Fertig, et al., Multiscale analysis of the effect of surface charge pattern on a nanopore's rectification and selectivity properties: from all-atom model to Poisson-Nernst-Planck, J. Chem. Phys., 150 (2019), 144703.
    [102] X. S. Wang, D. He, J. Wylie, H. Huang, Singular perturbation solutions of steady-state PoissonNernst-Planck systems, Phys. Rev. E, 89 (2014), 022722.
    [103] G. W. Wei, Q. Zheng, Z. Chen, K. Xia, Variational multiscale models for charge transport, SIAM Rev., 54 (2012), 699-754. doi: 10.1137/110845690
    [104] M. Zhang, Asymptotic expansions and numerical simulations of I-V relations via a steady-state Poisson-Nernst-Planck system, Rocky Mt. J. Math., 45 (2015), 1681-1708.
    [105] M. Zhang, Boundary layer effects on ionic flows via classical Poisson-Nernst-Planck systems, Comput. Math. Biophys., 6 (2018), 14-27. doi: 10.1515/cmb-2018-0002
    [106] J. Zhang, D. Acheampong, M. Zhang, Effects on ionic flows from finite ion sizes via PoissonNernst-Planck models with non-local excess chemical potentials, Comput. Math. Biophys., 5 (2017), 58-77. doi: 10.1515/mlbmb-2017-0005
    [107] Q. Zheng, G. W. Wei, Poisson-Boltzmann-Nernst-Planck model, J. Chem. Phys., 134 (2011), 194101.
    [108] Y. Zhou, R. MacKinnon, The occupancy of ions in the K+ selectivity filter: charge balance and coupling of ion binding to a protein conformational change underlie high conduction rates, J. Mol. Biol., 333 (2003), 965-975. doi: 10.1016/j.jmb.2003.09.022
    [109] M. Zhou, R. MacKinnon, A mutant KcsA K+ channel with altered conduction properties and selectivity filter ion distribution, J. Mol. Biol., 338 (2004), 839-846. doi: 10.1016/j.jmb.2004.03.020
  • This article has been cited by:

    1. Yiwei Wang, Lijun Zhang, Mingji Zhang, 2020, A special case study of boundary layer effects via Poisson-Nernst-Planck systems with permanent charges, 978-1-7281-7041-1, 1, 10.1109/ITNT49337.2020.9253312
    2. Mingji Zhang, Competition between Cations via Classical Poisson–Nernst–Planck Models with Nonzero but Small Permanent Charges, 2021, 11, 2077-0375, 236, 10.3390/membranes11040236
    3. Peter W. Bates, Zhenshu Wen, Mingji Zhang, Small Permanent Charge Effects on Individual Fluxes via Poisson–Nernst–Planck Models with Multiple Cations, 2021, 31, 0938-8974, 10.1007/s00332-021-09715-3
    4. Jianing Chen, Yiwei Wang, Lijun Zhang, Mingji Zhang, Mathematical analysis of Poisson–Nernst–Planck models with permanent charges and boundary layers: studies on individual fluxes, 2021, 34, 0951-7715, 3879, 10.1088/1361-6544/abf33a
    5. Jianing Chen, Mingji Zhang, Boundary layer effects on ionic flows via Poisson-Nernst-Planck systems with nonuniform ion sizes, 2022, 27, 1531-3492, 6197, 10.3934/dcdsb.2021312
    6. Lijun Zhang, Xiangshuo Liu, Chaohong Pan, Studies on reversal permanent charges and reversal potentials via classical Poisson-Nernst-Planck systems with boundary layers, 2022, 27, 1531-3492, 6631, 10.3934/dcdsb.2022013
    7. Jianing Chen, Mingji Zhang, Geometric Singular Perturbation Approach to Poisson-Nernst-Planck Systems with Local Hard-Sphere Potential: Studies on Zero-Current Ionic Flows with Boundary Layers, 2022, 21, 1575-5460, 10.1007/s12346-022-00672-0
    8. Zhenshu Wen, Peter W Bates, Mingji Zhang, Effects on I–V relations from small permanent charge and channel geometry via classical Poisson–Nernst–Planck equations with multiple cations, 2021, 34, 0951-7715, 4464, 10.1088/1361-6544/abfae8
    9. Mingji Zhang, Existence and local uniqueness of classical Poisson-Nernst-Planck systems with multi-component permanent charges and multiple cations, 2022, 0, 1937-1632, 0, 10.3934/dcdss.2022134
    10. Yanggeng Fu, Weishi Liu, Hamid Mofidi, Mingji Zhang, Finite Ion Size Effects on Ionic Flows via Poisson–Nernst–Planck Systems: Higher Order Contributions, 2022, 1040-7294, 10.1007/s10884-021-10114-1
    11. Yiwei Wang, Lijun Zhang, Mingji Zhang, STUDIES ON INDIVIDUAL FLUXES VIA POISSON-NERNST-PLANCK MODELS WITH SMALL PERMANENT CHARGES AND PARTIAL ELECTRONEUTRALITY CONDITIONS, 2022, 12, 2156-907X, 87, 10.11948/20210045
    12. Yanyu Bao, Jianing Chen, Lijun Zhang, Mingji Zhang, HIGHER ORDER EXPANSIONS IN FINITE ION SIZE VIA POISSON-NERNST-PLANCK SYSTEMS WITH BIKERMAN'S LOCAL HARD-SPHERE POTENTIAL, 2022, 12, 2156-907X, 907, 10.11948/20220001
    13. Xijun Deng, Yusheng Jia, Mingji Zhang, STUDIES ON CURRENT-VOLTAGE RELATIONS VIA POISSON-NERNST-PLANCK SYSTEMS WITH MULTIPLE CATIONS AND PERMANENT CHARGES, 2022, 12, 2156-907X, 932, 10.11948/20210003
    14. Yiwei Wang, Lijun Zhang, Mingji Zhang, Mathematical Analysis on Current–Voltage Relations via Classical Poisson–Nernst–Planck Systems with Nonzero Permanent Charges under Relaxed Electroneutrality Boundary Conditions, 2023, 13, 2077-0375, 131, 10.3390/membranes13020131
    15. Mingji Zhang, Qualitative properties of zero-current ionic flows via Poisson-Nernst-Planck systems with nonuniform ion sizes, 2022, 27, 1531-3492, 6989, 10.3934/dcdsb.2022030
    16. Yiwei Wang, Mingji Zhang, Finite ion size effects on I-V relations via Poisson-Nernst-Planck systems with two cations: A case study, 2024, 21, 1551-0018, 1899, 10.3934/mbe.2024084
    17. Hong Li, Zhantao Li, Chaohong Pan, Jie Song, Mingji Zhang, Cubic-like Features of I–V Relations via Classical Poisson–Nernst–Planck Systems Under Relaxed Electroneutrality Boundary Conditions, 2024, 13, 2075-1680, 790, 10.3390/axioms13110790
  • Reader Comments
  • © 2020 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(5123) PDF downloads(294) Cited by(17)

Figures and Tables

Figures(3)

Other Articles By Authors

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog