Geometrical effects on nonlinear electrodiffusion in cell physiology
aa r X i v : . [ q - b i o . S C ] M a y Geometrical effects on nonlinear electrodiffusion in cellphysiology
J. Cartailler , Z. Schuss ∗ , and D. Holcman † March 20, 2018
Abstract
We report here new electrical laws, derived from nonlinear electro-diffusion the-ory, about the effect of the local geometrical structure, such as curvature, on theelectrical properties of a cell. We adopt the Poisson-Nernst-Planck (PNP) equa-tions for charge concentration and electric potential as a model of electro-diffusion.In the case at hand, the entire boundary is impermeable to ions and the electricfield satisfies the compatibility condition of Poisson’s equation. We construct anasymptotic approximation for certain singular limits to the steady-state solution ina ball with an attached cusp-shaped funnel on its surface. As the number of chargeincreases, they concentrate at the end of cusp-shaped funnel. These results can beused in the design of nano-pipettes and help to understand the local voltage changesinside dendrites and axons with heterogenous local geometry.
Electro-diffusion is the process by which the motion of ions in solution is driven by twophysical forces: thermal motion, which is diffusion, and the electric field. The difficulty inthe mathematical description of this physical motion is due to the origin of the field, whichconsists of the contribution of mobile ions and of a possible external field. The dielectricmembrane also affects the field by image charges. So far only few electro-diffusion systemsare well understood: although the voltaic cell was invented more than 200 years ago,designing optimal configurations is still a challenge. On the other extreme, ionic flux andgating of voltage-channels [2] is now well explained by the modern Poisson-Nernst-Plancktheory of electro-diffusion [16], because at the nanometer scale, the cylindrical geometryapproximation of protein channels reduces the computation of the electric field and of ∗ Department of Mathematics, Tel-Aviv University, Tel-Aviv 69978, Israel. † Ecole Normale Sup´erieure, 46 rue d’Ulm 75005 Paris, France and Mathematical Institute, Uni-versity of Oxford, Andrew Wiles Building, Woodstock Rd, Oxford OX2 6GG, United Kingdom.Corresponding author email:[email protected] C and the surface S , that is, V ( C ) − V ( S ), increases, first linearly and then logarithmically, when the totalnumber of charges in the ball increases.Our aim here is to understand the effect of boundary curvature on an electrical cell,such as neuron. In particular, we explore the effect of boundary curvature on the chargeand field distribution at steady state. The curvature of membranes of dendrites andaxons of neurons have many local maxima that can modulate the channel’s local electricpotential [27]. In this article, we study the effects of local curvature on the distribution ofcharge in bounded domains with no electro-neutrality. The effect of non-electro-neutralitywas recently studied in [4] and a long-range electrostatic length, much longer than theDebye length was found. This effect is due to the combined effects of non-electro-neutralityand boundary, which lead to charge accumulation near the boundary.The cusp-shaped funnel geometry was studied in [10], however this paper presentsseveral crucial mathematical differences with [10], in particular, we are solving a nonlinearequation, while it was linear in [10]. Furthermore, the boundary condition at the end ofthe cusp-shaped funnel: while it is the Dirichlet condition in [10], it is the Neumanncondition here. This means that in [10] the absorption flux at the end of the funnel is2omputed, whereas here the stationary voltage and charge distribution are computed inthe absence of flux. We develop here new boundary layer analysis, different than theclassical matched asymptotics method [28, 29, 30]. The manuscript is organized as follow:first, we consider a bounded domain with an uncharged narrow cusp-shaped funnel onthe boundary, which is a singular geometrical effect. Second, we further study the case ofcharge distribution in a charged narrow cusp. The Poisson-Nernst-Planck system of equations in a domain Ω, whose dielectric boundary ∂ Ω is represented as the compatibility condition for Poisson’s equation, and its imper-meability to the passage of ions is represented as a no-flux boundary condition for theNernst-Planck equation. We assume that the total charge in Ω consists of N identicalpositive ions with initial particle density q ( x ) in Ω, their valence is z , and the total numberof particles is fixed, equal to Z Ω q ( x ) d x = N. (1)Thus the charge in Ω is Q = zeN, where e is the electronic charge. The charge density ρ ( x , t ) is the solution of the initialand boundary value problem for the Nernst-Planck equation D h ∆ ρ ( x , t ) + zekT ∇ ( ρ ( x , t ) ∇ φ ( x , t )) i = ∂ρ ( x , t ) ∂t for x ∈ Ω (2) D (cid:20) ∂ρ ( x , t ) ∂n + zekT ρ ( x , t ) ∂φ ( x , t ) ∂n (cid:21) = 0 for x ∈ ∂ Ω (3) ρ ( x ,
0) = q ( x ) for x ∈ Ω . (4)Here φ ( x , t ) is the electric potential in Ω and is the solution of the Neumann problem forthe Poisson equation ∆ φ ( x , t ) = − zeρ ( x , t ) εε for x ∈ Ω (5) ∂φ ( x , t ) ∂n = − σ ( x , t ) for x ∈ ∂ Ω , (6)where σ ( x , t ) is the surface charge density on the boundary ∂ Ω. In the steady state, σ ( x , t ) = Qεε | ∂ Ω | . (7)3 Steady solution in a ball with a cusp-shaped funnel
Local boundary curvature is a key geometrical feature that controls charge distributionin the domain. Specifically, we study the effect of a narrow funnel attached to a sphere.In various media, such as air (e.g., the lightning rod, [5]), the manifestation of this effectis observed in Lebesgue’s thorn, which is a an inverted cusp singularity of the boundary,for which the solution of Laplace’s equation blows-up inside the domain [5, p.304]. In thesteady state (2) gives the particle density ρ ( x ) = N exp (cid:26) − zeφ ( x ) kT (cid:27)Z Ω exp (cid:26) − zeφ ( x ) kT (cid:27) d x , (8)hence (5) gives Poisson equation∆ φ ( x ) = − zeN exp (cid:26) − zeφ ( x ) kT (cid:27) εε Z Ω exp (cid:26) − zeφ ( x ) kT (cid:27) d x . (9)and (6) gives the boundary condition ∂φ ( x ) ∂n = − Qεε | ∂ Ω | , (10)for | x | = R , which is the compatibility condition, obtained by integrating Poisson’s equa-tion (5) over Ω. Changing variables to u ( x ) = zeφ ( x ) kT , λ = ( ze ) Nεε kT , (11)Poisson’s equation (9) becomes∆ u ( x ) = − λ exp {− u ( x ) } Z Ω exp {− u ( x ) } d x (12)and the boundary condition (10) becomes ∂u ( x ) ∂n = − λ | ∂ Ω | for x ∈ ∂ Ω . (13)The translation ˜ u = u + ln (cid:18) λ/ Z Ω exp { v ( x ) } d x (cid:19) , converts (12) into − ∆˜ u ( x ) = exp {− ˜ u ( x ) } for x ∈ Ω (14) ∂ ˜ u ( x ) ∂n = − λ | ∂ Ω | for x ∈ ∂ Ω .
4e consider a dimensionless planar domain Ω with a cusp-shaped funnel formed by twobounding circles A and B of dimensionless radii 1 (see Fig.1(left)). The opening of thefunnel is ε ≪
1. We construct an asymptotic solution in this limit to the nonlinearboundary value problem (BVP) (14) by first mapping the domain Ω conformally with theM¨obius transformation of the two osculating circles A and B into concentric circles (seeFig.1(right)). To this end, we move the origin of the complex plane to the center of theosculating circle B and set w = w ( z ) = z − α − αz , (15)where α = − − √ ε + O ( ε ) . (16)The M¨obius transformation (15) maps the circle B (dashed blue) into itself and Ω ismapped onto the domain Ω w = w (Ω) in Figure 1(right). The straits in Figure 1(left) areFigure 1: Image Ω w = w (Ω) of the domain Ω ( A. ) under the conformal mapping (15). Theneck (left) is mapped onto the semi-annulus enclosed between the like-style arcs and the largedisk in Ω is mapped onto the small red disk. The short green segment AB (left) (of length ε ) is mapped onto the thick green segment AB (of length 2 √ ε + O ( ε )). The letters S and N designate the south and the north pole respectively. mapped onto the ring enclosed between the like-style arcs and the large disk is mappedonto the small red disk in Figure 1(right). The radius of the small disk and the elevationof its center above the real axis are O ( √ ε ). The short black segment AB of length ε inFigure 1(left) is mapped onto the segment AB of length 2 √ ε + O ( ε ) in Figure 1(right).This mapping (see [11]), transforms the PNP equations as well and thus leads to a new5on-linear effect. Setting u ( z ) = v ( w ) converts (12) to∆ w v ( w ) = − exp {− v ( w ) }| w ′ ( z ) | = − (4 ε + O (˜ ε / )) | w (1 − √ ˜ ε ) − O (˜ ε ) | exp {− v ( w ) } for w ∈ Ω w . (17)The boundary segment AB at the end of the cusp-shaped funnel in Figure 1(left) isdenoted ∂ Ω w,a . To determine the boundary conditions, we use the change of coordinates w = Re iθ = X + iY . At the end of the funnel, where R ≃
1, we get ∂u ( z ) ∂n z = − ∂v ( w ) ∂θ (cid:12)(cid:12)(cid:12)(cid:12) w = − ∂θ∂Y , (18)where ie iθ ∂θ∂Y = w ′ ( z ) = 1 − α (1 − αz ) . (19)For θ = π (for z = − ∂θ/∂Y = − / √ ε and the boundary condition at ∂ Ω w,a is ∂v ( w ) ∂n = − λ √ ε | ∂ Ω | for w ∈ ∂ Ω w,a . (20) Approximating the banana-shaped domain Ω w by a one-dimensional circular arc, we usea one-dimensional approximation of the solution in Ω w [12, 13]. This approximationassumes that there are no non-neutralized charges on the surface of the cusp (Fig.3A).The boundary condition for the approximate one-dimensional solution of (17) is zero atangle θ Lim = c √ ε , where c is a constant (see details in [12, 13]) and represents the solutioninside the disk in Figure 1(left), away from the cusp. Thus, (17) in the conformal imageΩ w becomes the boundary value problem v ′′ + 4 ε | e iθ − − e iθ √ ε | exp (cid:8) − v ( e iθ ) (cid:9) = 0 (21) v ′ ( c √ ε ) = 0 (22) v ′ ( π ) = − λ √ ε | ∂ Ω | . The solution of (21) is shown in Figure 3B-C in the two domains, Ω (panel A ) and itsimage Ω w (panel B ). 6 x y x -3 . . . y x y A B
Figure 2:
Influence of the cusp on the field lines (orthogonal to the level lines) . Thefield line inside the original domain Ω (A) and its image domain Ω w ( B), computed numericallyfrom equation (14). The blue lines originate from the bulk, while the orange starts in the cusp.The domain Ω w is subdivided into three regions: the region Ω w inside the funnel, the region Ω w connecting the end of the funnel to the bulk Ω w . Our goal is now to estimate the difference of potentials between the north pole N andthe end of the funnel C ,˜∆ u = u ( N ) − u ( C ) = v ( c √ ε ) − v ( π ) . (23)To construct an asymptotic approximation to the solution of (22) in the limits ε → λ → ∞ , we first construct the outer-solution in the form of a series in powers of ε , whichis an approximation valid away from the boundary. In the limit of small ε , the first termin the series vanishes, exponential terms drop out, and the second order term is y outer ( θ ) = M θ + M ′ , (24)where M and M ′ are yet undetermined constants. The outer solution cannot satisfy allboundary conditions, so a boundary layer correction is needed at the reflecting boundaryat θ = c √ ε . Thus, we set θ = √ εξ and expand ε | e iθ − − e iθ √ ε | = 1(1 + ξ ) + O ( √ ε ) . Writing the boundary layer solution as y bl ( θ ) = Y ( ξ ), we obtain to leading order theboundary layer equation Y ′′ ( ξ ) + 4(1 + ξ ) exp {− Y ( ξ ) } = 0 , (25)7ith Y ′ ( c ) = 0. The solution is decaying for large ξ and develops a singularity at finite ξ .However, a Taylor expansion near ξ = 0, Y ( ξ ) = A + B ξ + B ξ + . . . , (26)gives in (25) B = − e − A . (27)In general, the coefficients satisfy B k = O ( e − A ), for A ≫
1. For small ξ , we obtain theapproximate solution of (25) by considering the leading term in a regular expansion ofthe solution in powers of ξ . The equation for the leading term is Y ′′ ( ξ ) + 4 e − A (1 + ξ ) = 0 (28)and the solution is defined up to an additive constant. Setting Y appr (0) = 0, which doesnot affect the potential difference, we find that Y appr ( ξ ) = − ξe − A arctan ξ. (29)It follows that the boundary layer solution at c √ ε is y bl ( θ ) = A − θ √ ε e − A arctan θ √ ε . (30)The boundary layer near π is needed, because A → ∞ as ε → b > d dθ v ( θ ) + be − v ( θ ) = 0 , dv (0) dθ = v (0) = 0 . The solution is v b ( θ ) = ln cos b θ. (31)Putting the outer and boundary layer solutions together gives the uniform asymptoticapproximation y unif ( θ ) = A − θ √ ε e − A arctan θ √ ε + ln cos b θ, (32)where the parameters A and b are yet undetermined constants. The condition at c √ ε = o (1) for ε ≪ y ′ unif (0) = 0 . θ = π gives that y ′ unif ( π ) = − πe − A √ ε − b tan b π = − λ √ ε | ∂ Ω | . The compatibility condition for (14), λ = Z Ω exp {− ˜ u ( x ) } dS x , (33)gives in Ω w that λ = Z Ω w exp {− ˜ v ( w ) } dw | φ ′ ( φ − ( w )) | = 8 √ ε π Z c √ ε exp {− v ( θ ) }| e iθ (1 − √ ε ) − | dθ. (34)Using the uniform approximation (32) in the compatibility condition (34), we obtain thesecond condition λ = 8 √ εe − A π Z c √ ε b θ exp (cid:26) e − A θ √ ε arctan θ √ ε (cid:27) | e iθ (1 − √ ε ) − | dθ ≈ e − A ε π/ √ ε Z b √ εξ exp (cid:8) e − A ξ arctan ξ (cid:9) | ξ | dξ, (35)where we used the change of variable θ = √ εξ . Integrating by parts, we get for ε ≪ λ ∼ e − A ε b √ ε tan b π exp (cid:26) e − A π √ ε π (cid:27)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (cid:18) π √ ε (cid:19) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) − π/ √ ε Z b √ ε tan b θ Ψ( θ ) dθ , (36)where Ψ( ξ ) = ddξ exp (cid:8) e − A ξ arctan ξ (cid:9) | ξ | . (37)Thus, it remains to solve the asymptotic equation λ ∼ e − A ε / (cid:20) bπ tan πb (cid:26) π e − A √ ε (cid:27) + O (cid:18) ln (cid:12)(cid:12)(cid:12)(cid:12) cos πb (cid:12)(cid:12)(cid:12)(cid:12)(cid:19)(cid:21) . (38)9or A and b in the limit ε →
0. We consider the limiting case where e − A √ ε = O (1) = C for λ → ∞ , (39)for which condition (33) can be simplified and gives to leading order b tan πb λ √ ε | ∂ Ω | , (40)that is, for λ √ ε ≪ b ≈ − π | ∂ Ω | λ √ ε , tan b π ∼ λ √ ε | ∂ Ω | . With condition (38), we get λ ≈ e − A ε / (cid:20) π λ √ ε | ∂ Ω | exp (cid:26) π e − A √ ε (cid:27) + O (cid:18) ln (cid:12)(cid:12)(cid:12)(cid:12) cos πb (cid:12)(cid:12)(cid:12)(cid:12)(cid:19)(cid:21) , (41)To leading order in large C , we obtain π | ∂ Ω | ε / = C exp (cid:8) Cπ (cid:9) . (42)The solution is expressed in terms of the Lambert-W function, Cπ = W (cid:18) π | ∂ Ω | ε / (cid:19) , (43)and for small ε , using the asymptotics of the Lambert function, Cπ = ln π | ∂ Ω | ε / − ln (cid:20) ln π | ∂ Ω | ε / (cid:21) + o (1) . (44)Finally, e − A √ ε = C ∼ π ln π | ∂ Ω | ε / ,A = ln 1 √ ε − ln (cid:20) π ln π | ∂ Ω | ε / (cid:21) → ∞ as ε → . (45)It follows that a uniform asymptotic approximation (32) in the limits λ → ∞ ε → y unif ( θ ) = ln 1 √ ε − ln (cid:20) π ln π | ∂ Ω | ε / (cid:21) (46) − θ π ln π | ∂ Ω | ε / arctan θ √ ε + ln " cos − π | ∂ Ω | λ √ ε θ . ε and λ in Figure 3against the numerical solution of (21), with the boundary conditions v ′ ( c √ ε ) = v ′ (0) = 0.The numerical solutions are computed with the software COMSOL, based on an adaptivemesh refinement and a relative tolerance of 10 − , that we validated on known analyticalresults of steady state PNP equations in a disk [4]. We find that the asymptotic expansionis particularly good in the limit ε → λ → ∞ (Fig.3A-D). However, for λ = O (1)the log-term approximation in (46) is non-monotonic in θ . Finally, to further validate the B A Reduced PNP:Approx.
CD E V nu m - Figure 3:
The asymptotic solution y unif ( θ ) of (32) (blue dashed lines) is compared to thenumerical solution of (21) (red line). The four panels A-B-C-D are obtained for different pairsof parameters ( λ, ε ). E. shows a 3D plots the difference between the asymptotic solution y unif (eq. 32) and numerical results V num (eq. 21), averaged over the domain Ω w .uniform asymptotic expansion, we compared the numerical solutions of the full equation(23) in the initial domain Ω with the reduced PNP equation (14) with zero Neumannboundary conditions, except at the end of the funnel for the mapped domain Ω w . Theresult is shown in Figure 3A-D, showing good agreement between the one-dimensionalPNP approximation in Ω w and the numerical solution of the full equation. We estimatednumerically the difference between the asymptotic solution y unif (32) and the numericalestimation V num (21), averaged over the domain Ω w , for 10 ≤ λ ≤ · and 5 · − ≤ ε ≤ − . The difference is almost constant in the range [0 . , . u ( N ) − u ( S ) has a maximumwith respect to λ , where u ( N ) and u ( S ) are the values of the potential at the north poleand at the end of the funnel, respectively. 11 .2 The voltage drop between the end of the funnel and thecenter of the ball We can now use (32) to compute the potential drop in (23). It is given by˜∆ SC u = u ( S ) − u ( C ) = − v ( c √ ε ) + v ( π )= − ln π | ∂ Ω | ε / + 2 ln 2 | ∂ Ω | λε / = ln 2 | ∂ Ω |√ επ λ . (47)The potential difference ˜∆ SC u with respect to λ is shown in Figure 3F (red line).Next, we compare the potential drop (23) with the one between the center and thenorth pole. Numerical solution of the PNP equations shows that the voltage and chargedistribution in a disk with a funnel do not differ from the ones in a disk in the uppersphere (Fig.2). This result is compared next to the difference between the north pole andthe center evaluated from the exact analytical expression derived for a disk.The expression for the voltage in the two-dimensional disk of radius R is given by (see[4]) u Dλ ( x ) = ln (cid:20) − λ D π + λ D (cid:16) rR (cid:17) (cid:21) , where λ D is a parameter. We calibrate λ D so that the solutions of the PNP equationsin a disk with a funnel have the same total charge as a disk. The Neumann boundaryconditions for the disk and the funnel are, respectively, ∂u ( x ) ∂n = − λ D πR , ∂u ( x ) ∂n = − λ | ∂ Ω | . The calibration is λ D = λ πR | ∂ Ω | . (48)We compare in Figure 3D the two-dimensional numerical solution of the PNP equation(14) in the domain Ω (blue line), with the analytical solution (48) in a disk with no cusp(dashed red). The numerical solution of the PNP equation (14) is plotted along the mainaxis 0 y in the interval [0 , y ] (where the point y is defined by the condition ∇ u ( y ) = 0).In the range [ y , y cusp ], where y cusp is the coordinate of the cusp, we compare the solutionof (14) with the uniform solution y unif of (32) in the funnel (dashed green). We concludethat in the cusp, the two-dimensional approximation in a disk is in good agreement withthe numerical solution of equation (14), confirming that the solution in the bulky headdoes not influence the one in the cusp (as already shown in Fig. 2). This result alsoconfirms the validity of the analytical formula to predict the large λ asymptotics.For a disk of radius R , the potential drop is given by˜∆ NC u = u ( N ) − u ( C ) = ln (cid:18) π π + λ D (cid:19) = − λ − (cid:18) R | ∂ ˜Ω | (cid:19) + O ( 1 λ ) (49)12see section 3.2). The potential drop ˜∆ NC u is shown in Figure 3E (blue line). The twodifferences of potential ¯∆ SC u (49) and ˜∆ NC u (47) have the same logarithmic behaviorln 1 /λ for λ ≫ u ( N ) − u ( S ) = O (1). A numerical solution in two-dimensionsshows that u ( N ) − u ( S ) may converge to zero as λ increases (Fig. 3F), thus having alocal maximum for small values of λ . This maximum cannot be analyzed by the uniformexpression (32), because it appears outside the domain of validity of (32). This resultis in agreement with the two-dimensional numerical solution of (14) for the differencebetween u ( N ) (potential at the north pole) and u ( S ) (potential at the end of the funnel)(Fig. 3F). The potential drop calculated above is non-dimensionalized by the radius ofcurvatures R f at the right and left of the funnel, ε = ˜ εR f , where ˜ ε is the length of the absorbing arc AB . The non-dimensionalized volume andboundary measure are, respectively, | Ω | = | ˜Ω | R , | ∂ Ω | = | ∂ ˜Ω | R .
In dimensional units (47) gives the potential drop in the dimensional disk with a funnelas ˜∆ SC u = u ( S ) − u ( C ) = ln 2 | ∂ ˜Ω |√ ˜ επ R / f λ . (50)We conclude in the limit of λ ≫
1, ˜ ε → S and the north pole N in the domain is obtained by adding (49) and (50) andwe get˜∆ SN u = u ( S ) − u ( C ) + u ( C ) − u ( N ) = ln | ∂ ˜Ω |√ ˜ επ R / f ! + 2 ln (cid:18) R | ∂ ˜Ω | (cid:19) + O ( 1 λ ) . (51)We recall that R is the radius of the entire ball, while R f is the radius of curvature of thefunnel. Due to the Neumann boundary condition (6) on the lateral part of the funnel,(17) in thetransformed domain cannot be reduced to one dimension. Thus we derive a different one-dimensional approximation for the mapped PNP equations in the banana-shaped domain13 A y u ( y )
2D numericsDisk analyticalReduced PNP D y y Figure 4:
Comparison of the numerical solutions of the full and reduced PNP equations(14) with zero Neumann boundary conditions, except at the end of the funnel. A. Schematicrepresentation of the domain Ω with an uncharged cusp (blue). The letters N , S , and C referto the north pole, the funnel tip, and the center of mass respectively. B-C
Numerical solutionsof (14) (solid) and the solution of (57) in the funnel (dashed) in the mapped domain Ω w . Thesolution have been obtained for ε = 0 . D. Comparison of (14) (blue) with the numericalsolution (21) inside the funnel (dashed green) and (48) in the bulk (dashed red). E. Solution u ( S ) − u ( C ) (dashed blue) obtained numerically from (47) and compared to the logarithmicfunction − λ ) (greed dotted). F. Two-dimensional numerical solutions of the difference | u ( N ) − u ( C ) | vs λ . The inset in panel F. is a blowup showing a maximum for small λ . w by averaging over the radius r . Rewriting (17) in polar coordinates w = re iθ , weobtain1 r ∂∂r (cid:18) r ∂v ( w ) ∂r (cid:19) + 1 r ∂ v ( w ) ∂θ = − (4 ε + O ( ε / )) exp {− v ( w ) }| re iθ (1 − √ ε ) − O ( ε ) | for w ∈ Ω w . (52)In the section Ω w ∩ { − √ ε < r < } , the boundary conditions are ∂v ( r, θ ) ∂r (cid:12)(cid:12)(cid:12)(cid:12) r =1 = − λ √ ε | ∂ Ω | (cos θ − , for θ ∈ [ c √ ε, π ] (53) ∂v ( r, θ ) ∂r (cid:12)(cid:12)(cid:12)(cid:12) r =1 −√ ε = 0 , for θ ∈ [ c √ ε, π ] ∂v ( r, θ ) ∂θ (cid:12)(cid:12)(cid:12)(cid:12) θ = π = − λ √ ε | ∂ Ω | ,∂v ( r, θ ) ∂θ (cid:12)(cid:12)(cid:12)(cid:12) θ = c √ ε = 0 . Taylor’s expansion of v in the section gives v ( r, θ ) = v ( θ ) + ( r − v ( θ ) + O (( r − ) , (54)and because | r − | = O ( √ ε ), we obtain the approximation,exp {− v ( w ) } = exp {− v ( θ ) } (cid:0) − √ εv ( θ ) + O ( ε ) (cid:1) . Multiplying (52) by r and integrating over the radius, we get (cid:20) r ∂v ( r, θ ) ∂r (cid:21) r =11 −√ ε + ∂ ∂θ Z −√ ε v ( r, θ ) dr = − Z −√ ε (4 r ε + O ( ε / )) e − v ( r, θ ) | re iθ (1 − √ ε ) − O ( ε ) | dr. (55)The boundary conditions (53) give, to leading order in √ ε , that − λ √ ε | ∂ Ω | (cos θ −
1) + √ ε ∂ v ( θ ) ∂θ = (56) − Z −√ ε (4 r ε + O ( ε / )) | re iθ (1 − √ ε ) − O ( ε ) | e − v ( θ ) (cid:0) − √ εv ( θ ) + O ( ε ) (cid:1) dr. that is, the BVP (52) in the section becomes the ODE (with respect to θ ), v ′′ ( θ ) = − (4 ε + O ( ε / )) | e iθ (1 − √ ε ) − O ( ε ) | exp {− v ( θ ) } − λ | ∂ Ω | (1 − cos θ ) , (57) v ′ ( θ ) | θ = π = − λ √ ε | ∂ Ω | ,v ′ ( θ ) | θ = c √ ε = 0 . ω and Ω is shown in Figure 6. Equation (57)is obtained by averaging over the radial direction and its solution seems to be a goodapproximation to (52) only for small λ . A different approach for large λ is discussed inthe next section.A regular expansion for λ ≪ v ( θ ) = w ( θ ) + λw ( θ ) + o ( λ ) , (58)gives in (57) that w = O ( ε ) and w is the solution of the BVP w ′′ ( θ ) = − (4 ε + O ( ε / )) | e iθ (1 − √ ε ) − O ( ε ) | − | ∂ Ω | (1 − cos θ ) , (59) w ′ ( θ ) | θ = π = −√ ε | ∂ Ω | , (60) w ′ ( θ ) | θ = c √ ε = 0 . Direct integration with respect to θ gives w ( θ ) = − θε √ ε arctan θ √ ε + 1 | ∂ Ω | ln sin θ Aθ + B. (61)Equation (60) gives A as A = πε / − π − √ ε | ∂ Ω | . (62)The zero Neumann boundary condition cannot be satisfied and a boundary layer appears,leading to the local expansion v ( θ ) = λw ( θ ) + o ( λ ) . (63)It follows that for λ ≪
1, the solution increases with λ . It is shown below that it decreasesfor λ ≫
1, demonstrating that there is at least one maximum in the variable λ . In the limit of λ ≫ ε →
0, the asymptotic expansion of the potential found above for acharged disk with a funnel is no longer valid. Some insight can be gained by observing thefield lines in the domain Ω w , described in Figure 2B. These lines are parallel to the radiusvector, except in a small region near the funnel. Two sections can be distinguished, A = { ( r, θ ) ∈ Ω w : | θ − √ ε | > π, | r − | ≤ √ ε } (64) B = { w = (1 − √ ε ) e iθ : | θ − π | ≤ √ ε } . The two sections A and B are illustrated in Figure 5A. Note that the boundary of section B contains a circular arc (marked magenta). Next, the approximate solutions u A ( r, θ )and u B ( θ ) of (52) in the two sections are constructed and used to construct a uniformapproximation u unif in Ω w (Fig. 5B). 16igure 5: Decomposition of the banana-shaped domain Ω w into two subregions regions A and B . A. Representation of the two subregions A (blue) and B (magenta) of Ω w . B. Solutionsof (83) (dashed blue), (90) (red dots), and the uniform approximation u unif of (32) (green) for r = 1 − √ ε . u A ( r, θ ) in section A The boundary conditions (53) for the potential equation (52) indicate that the radialderivative is O ( λ √ ε ) → ∞ Thus the angular derivatives are negligible relative to theradial ones. It follows in a regular expansion of the solution that the θ derivatives canbe neglected relative to the r derivative and the equation is then solved along the rays θ = const = θ for r ∈ [1 − √ ε, λ √ ε , u ′′ A ( r, θ ) + 1 r u ′ A ( r, θ ) = − ε exp( − u A ) | re iθ (1 − √ ε ) − | for r ∈ [1 − √ ε,
1] (65) u ′ A ( r, θ ) | r =1 −√ ε = 0 u ′ A ( r, θ ) | r =1 = − λ √ ε | ∂ Ω | (1 − cos θ ) . For ε ≪
1, we get | re iθ (1 − √ ε ) − | = | e iθ − | + O ( √ ε ). Setting h ( θ ) = 4 ε | e iθ − | , (66)and v A,θ ( r ) = − u A ( r, θ ) + ln h ( θ ) , (67)17e get v ′′ A,θ ( r ) + 1 r v ′ A,θ ( r ) = exp( v A,θ ) (68) v ′ A,θ ( r ) (cid:12)(cid:12) r =1 −√ ε =0 v ′ A,θ ( r ) (cid:12)(cid:12) r =1 = λ √ ε | ∂ Ω | (1 − cos θ ) . The general solution of (68) is given by [4] v A,θ ( r ) = ln C r − ln cos C r − C ) , (69)where the constants C and C are determined from the boundary conditions (68). Using v ′ A,θ ( r ) = C r tan C r − C ) − r , (70)we find the constant C from (70) and from the boundary condition (68) at the point r = 1 − √ ε , getting C = − (cid:18) C arctan 2 C + √ ε (cid:19) + O ( ε ) . (71)This gives in (70) at r = 1 the transcendental equation for C , C tan − C C λ √ ε | ∂ Ω | (1 − cos θ ) + 2 , (72)hence lim λ →∞ − C C π . (73)Now, it follows from (71) that − C C C + C √ ε. (74)Note that lim λ →∞ C = 0, because otherwise we would get the asymptotic expansion − C C π C √ ε −
1) + O ( C ) , (75)which leads to C tan − C C − √ ε + O ( C ) (76)18nd contradicts the condition (72) in the limit λ → ∞ .Then (73) and (74) would imply that C √ ε O (1) (77)and (77) would give C ≫
1, so that the arctan term in (74) drops out, and we would beleft with − C C ∼ C √ ε, (78)hence C ∼ −√ ε. (79)Expanding the left hand side of (72), using (73) and (78), we obtain thattan C √ ε − C √ ε − π + O (cid:18) C √ ε − π/ (cid:19) . (80)Together with (80), the solution of (72) is C ∼ λπ √ ε | ∂ Ω | (1 − cos θ ) + λε . (81)With the values of C and C computed in (71) and (81), the solution v A,θ of (69) isgiven by v A,θ ( r ) = ln ε r (cid:18) λπ | ∂ Ω | (1 − cos θ ) + λε (cid:19) (82) − ln cos λπ √ ε [ln r + √ ε ]2 | ∂ Ω | (1 − cos θ ) + λε . Finally, using (67) and (82), we obtain for ( r, θ ) ∈ A , u A ( r, θ ) = − ln | e iθ − | r (cid:18) λπ | ∂ Ω | (1 − cos θ ) + λε (cid:19) (83)+ ln cos λπ √ ε [ln r + √ ε ]2 | ∂ Ω | (1 − cos θ ) + λε . The asymptotic solution u A is plotted in Figure 5B (blue dashed line). Comparison withnumerical solutions for various values of λ and ε is shown in Figure 6 below.19 .3 The asymptotics of u B in section B The asymptotic solution u A ( r, θ ) in section A cannot satisfy the boundary condition (53)at θ = π . Indeed, (83) gives ∂u A ( r, θ ) /∂θ | θ = π = 0, while the boundary condition (57) is ∂v/∂θ | θ = π = − λ √ ε/ | ∂ Ω | , so a boundary layer correction is needed.The boundary layer u B ( θ ) is an asymptotic solution of (52) in section B , where the θ derivatives dominate the radial ones. The right-hand-side of (52) can be simplified for ε ≪ r = 1 − √ ε the approximation − ε | re iθ (1 − √ ε ) − | ∼ − ε r and θ . With this simplification in (52), we rewrite u B ( θ ) as u B ( θ ) = ˜ u B ( η ) + C , (85)where C is an additive constant and ˜ u B is a function of η = θ − ( π − √ ε ) and solves theBVP ∂ ˜ u B ( η ) ∂η = − exp {− ˜ u B ( η ) } (86)˜ u ′ B ( η ) | η = √ ε = − λ √ ε | ∂ Ω | ˜ u ′ B ( η ) | η =0 =0 . The solution of (86) (see [4]) is˜ u B ( η ) = ln cos r λ I λ η, (87)where I λ is the solution of the transcendental equation I λ = 2 | ∂ Ω | λε tan r λε I λ . (88)We obtain to leading order for λ ≫ I λ = 2 λεπ (cid:18) | ∂ Ω | λε (cid:19) + O (cid:18) λε (cid:19) . (89)It follows from (89), (87), and (85) that for θ ∈ B , the asymptotic solution is u B ( θ ) = ln cos π r ( θ − ( π − √ ε )) ε (cid:18) − | ∂ Ω | λε (cid:19) + C (90)(see (85)). It is shown in Figure 5B (red dots).20 .4 A uniform approximation of u ( r, θ ) in Ω w A uniform asymptotic approximation u unif ( r, θ ) of the voltage u ( r, θ ) in the entire mappeddomain Ω w can be now constructed by matching the the leading term u A ( r, θ ), given in(83) in section A , with that of u B ( θ ), given in (90) in section B .These approximations agree at θ = π − √ ε , so we obtain that C = u A (1 − √ ε, π − √ ε ) . (91)Thus u unif ( r, θ ) = ( u A ( r, θ ) for θ ∈ [0 , π − √ ε ] u B ( θ ) for θ ∈ [ π − √ ε, π ] . (92)The numerical solution of (14) in Ω w and the approximation u unif ( r, θ ) of (92) are shownFig. 5C. Ω w The potential drop ˜∆ funnel u between the center of mass C and the tip of the funnel S , is∆ funnel u = u ( C ) − u ( S ) . (93)Due to the axial symmetry of the domain Ω, the center of mass C is at r = 1 − √ ε , hence(92) gives u ( S ) = u (1 − √ ε, π ) and u ( C ) = u (1 − √ ε, c √ ε ) . (94)Recall that the constant c depends on the domain geometry only, and is defined by theconformal mapping w (see relation (15)). The potential drop ˜∆ Cusp u in the funnel canbe decomposed as the sum of difference of potential between the two sections, A and B .First, the approximations are˜∆ u A = u A (1 − √ ε, π ) − u A (1 − √ ε, c √ ε ) . (95)and ˜∆ u B = u B ( π ) − u B ( π − √ ε ) , (96)so that ˜∆ funnel u ∼ ˜∆ u A + ˜∆ u B . (97)Using (83) in A , we get u A (1 − √ ε, θ ) = − ln | e iθ − | − √ ε ) (cid:18) λπ | ∂ Ω | (1 − cos θ ) + λε (cid:19) − ln cos λπ √ ε (ln(1 − √ ε ) + √ ε )2 | ∂ Ω | (1 − cos θ ) + λε . (98)21or ε ≪
1, we get from (98) that − ln cos λπ √ ε (ln(1 − √ ε ) + √ ε )2 | ∂ Ω | (1 − cos θ ) + λε = O ( ε ) . (99)Hence, using (99) in (98), we get u A (1 − √ ε, θ ) = − ln | e iθ − | − √ ε ) (cid:18) λπ | ∂ Ω | (1 − cos( θ )) + λε (cid:19) + O ( ε ) . (100)The approximate solution u A ( S ) at the tip of the funnel S (south pole at θ = π ) is(100) u A ( S ) = − ln 2 λ π (4 | ∂ Ω | + λε ) + 2 ln(1 − √ ε ) + O ( ε ) . (101)At the center C , where θ = c √ ε , equation (98) gives for ε ≪ θ -dependent termsin (100) as | e iθ − | = c ε + O ( ε ) , (102)and 2 | ∂ Ω | (1 − cos c √ ε ) + λε = ε ( | ∂ Ω | c + λ ) + O ( ε ) . (103)Using (102) and (103), the expression (98) reduces to u A ( C ) = − ln c (cid:18) λπ | ∂ Ω | c + λ (cid:19) + 2 ln(1 − √ ε ) + O ( ε ) . (104)For λ ≫
1, (104) becomes u A ( C ) = − ln π c − √ ε ) + O (cid:18) ε, λ (cid:19) red. (105)Finally, the approximate potential difference ˜∆ u A in (95), is the difference between (105)and (101), ˜∆ u A = − ln 2 λ π (4 | ∂ Ω | + λε ) + ln π c O (cid:18) ε, λ (cid:19) . (106)For λ ≫ u A ∼ − ln 2 c ε , (107)22hich is independent of λ . (90) shows that the approximate potential in section B is u B ( π − √ ε ) = C (108)and u B ( π ) = ln sin π | ∂ Ω | λε + C . (109)Using (108) and (109) in (96), we obtain˜∆ u B = ln sin π | ∂ Ω | λε . (110)For λ ≫
1, (110) shows that ˜∆ u B is˜∆ u B = − λ + 2 ln | ∂ Ω | πε + O (cid:18) λ (cid:19) . (111)Finally, using (106), (110) and (97), we find that the potential drop is˜∆ u = ln sin π | ∂ Ω | λε − ln 2 λ π (4 | ∂ Ω | + λε ) + ln π c O (cid:18) ε, λ (cid:19) . (112)Again, using (107), (111) and (97) for λ ≫ u ∼ − ln λ + 2 ln πc | ∂ Ω | O (cid:18) λ (cid:19) . (113)Equation (110) shows that for λ ≫
1, the potential drop in the funnel domain occursmostly in the region B . The expression (112) is plotted in Figure 5A-D (red) and comparedto ln λ + const (green) and to a two-dimensional numerical solution. The good agreementconfirms the validity of the asymptotic expansion and thus confirming the new asymptoticformulas derived here. We conclude with the general formula for a dimensional cusp-shaped funnel where | ∂ Ω | = | ∂ ˜Ω | R c and R c is the radius of curvature at the cusp˜∆ u ∼ − ln λ + 2 ln πc | ∂ ˜Ω | R c + O (cid:18) λ (cid:19) . (114)23 F B C
2D numerics2l og ( ) +C ste u uni f E
2D numerics u uni (cid:0) - u ---------- Figure 6:
Comparison of numerical solution of (14) in the plane with the approximations u unif ( x ) in (92). A. Schematic representation of the domain Ω with a charged funnel (red). Theletters N , S , and C refer to the north pole, the funnel tip, and the center of mass, respectively. B-C
Numerical solutions of (14) (solid) and the solution of (92) in the funnel (dashed) in themapped domain Ω w for several values of λ and for ε = 0 . D. Comparison of (14) (blue) withanalytical solutions (32) inside the funnel (dashed green) and (48) in the bulk (dashed red). E. Solution u ( S ) − u ( C ) obtained numerically (dashed blue) from (47) and analytically from (32)(red), compared to the logarithmic function − λ + const (green dots). F. Two-dimensionalnumerical solutions of the difference | u ( N ) − u ( C ) | vs λ compared to the analytical solutions(112) (red). The inset in panel F. is a magnification showing a maximum for small λ . .6 Expansion of the potential drop between N and S To expand the potential difference u ( N ) − u ( S ) between the funnel tip S and the northpole N of Ω, we first use the results (114) computed above, to expand the difference u ( C ) − u ( S ), and then subtract (114) and (49). The the terms 2 ln( λ ) drop out and wehave u ( N ) − u ( S ) = 2 ln 4 | ∂ Ω | R − πc | ∂ Ω | R c + O (cid:18) λ (cid:19) , (115)where R is the distance between the north pole N and the center of mass C and R c is theradius of curvature at the cusp. We obtain to leading order u ( N ) − u ( S ) ∼ − πc R R c , (116)which is a constant that depend only on the center of mass C . A C (cid:1) harged (cid:2) us p Uncharged Cusp N o r m a li ze d B y Figure 7:
Normalized charge distribution ρ ( y ) in charged and uncharged funnel domains. A. ρ ( y ) is computed numerically from (8) with ∂u/∂n = 0 at the funnel boundaries ( λ = 1 (blue), λ = 10 (red), λ = 1500 (green), and λ = 1000 (dashed magenta)). B. Representation of Ω andthe funnel boundary conditions. Left: uncharged funnel domain ∂u/∂n = 0 (blue), and Right:charged funnel domain ∂u/∂n = − λ/ | Ω | (red). C. ρ ( y ) in a charged funnel domain. The samecolor code is used as in panel A. We have derived here new electrostatic laws in non-neutral confined electrolytes fromnonlinear electro-diffusion theory (PNP equations). The effect of local geometrical struc-ture, such as the local curvature of the boundary emerges from the asymptotic solution of25he model. The PNP equations describe the charge concentration and electric potential.The new electrical laws are derived in the context of non-electro-neutrality, where we usea single ionic species. The approximation of the steady-state solution in a ball with anattached cusp-shaped funnel on its surface is new and the construction of the asymptoticexpansion uses a new boundary layer analysis.Using asymptotic and numerical solution of the PNP equation, we found here that fora sufficiently high number of charges, the charge concentration peaks at the end of thefunnel in a charged funnel boundary domain, but this is not the case for an unchargedfunnel domain (Fig.7A-C). This effect is clearly the result of the cusp-shaped geometry.The present analysis reveals that the curvature affects the membrane potential. We alsofind that the voltage increases logarithmically in the total number of excess charges N ,which is valid for uncharged (47) and charged (52) cusp-shaped funnel on the boundary.We studied here the voltage changes and electro-diffusion under an excess of positiveions. The voltage difference in the limit λ → ∞ is probably attenuated in a mixed ionicsolution, but the electro-neutrality remains broken. Cytoplasmic ions are characterizedby the following concentrations Na + = 148ml, K + = 10ml and Cl − = 4ml. There is a clearunbalance toward positive charges, however there are probably molecules of various sizeswith negative charges to re-balance the charges. However, the motility of these proteinsshould be driven by a diffusion coefficient smaller than the one of the ions. This differenceof mobility is certainly a key feature in maintaining non-electro-neutrality and then tuningthe value of λ . However, in a system containing an excess of positive and a small amountof negative charges, we show in Appendix that the limit of the PNP equations in thebulk, when the number of negative charges tends to zero, can be obtained by a regularexpansion of the solution. This result shows that the small amount of negative chargesdoes not perturb much the distribution of positive ones. red Finally, note that we did notconsider here nanometer structures, such as ionic channels, where a negative ionic chargecan affect the motion of the other ions in the channel pore.We conclude that local geometrical properties, such as curvature, can modulate thelocal voltage in biological cellular electrolytes when electro-neutrality is violated. Thisresult generalizes the case of a ball, where the distribution of charges accumulates on thesurface as the total charge increases [4]. Following a non uniform boundary curvature, weexpect that charges will be non-uniformly distributed, leading to a difference of potentialacross the membrane with charges on its surface. Since, this difference of potential playsa key role in information processing at synapses, we conclude that the spine geometry,in particular its curvature may impact the coding or decoding of voltage through current[27].This effect may as well influence the propagation and genesis of local depolarization[20, 19, 14]. More realistic funnels, with two different curvature radii can be incorporatedto the formalism presented by modifying the parameter α (15) as shown in [12]. Theformalism presented in this paper can be applied beyond physiology, in particular in thedesign of nanopipettes with an optimal shape [18, 12] by modulating α (15) or with apatterned surface [26] by changing the surface charge density via λ in region A (65).26 Appendix red
Regular expansion of the PNP solution when there are an excessof positive and a small number of negative charges
We show that for the concentrations of cations and anions found in literature [31], theleading order solution of the electrical potential in the bulk can be obtained by consideringpositive charges only. We assume that the charge of an electrolyte confined in ˜Ω consistsof identical N p positive and N m negative ions with density q p ( x ) and q m ( x ) such as N i = Z ˜Ω q i (˜ x ) d ˜ x , for i ∈ { p , m } , (117)where p and m are positive and negative species respectively. The total charges in ˜Ω isthe sum Q = e ( N p − N m ) . (118)The associated charge densities ρ p ( x , t ) and ρ m ( x , t ) satisfy the boundary value problemfor the Nernst-Planck equation D i h ∆ ρ i (˜ x , t ) + z i ekT ∇ ( ρ i (˜ x , t ) ∇ φ (˜ x , t )) i = ∂ρ i (˜ x , t ) ∂t for ˜ x ∈ ˜Ω (119) D i (cid:20) ∂ρ i (˜ x , t ) ∂n + z i ekT ρ i (˜ x , t ) ∂φ (˜ x , t ) ∂n (cid:21) = 0 for ˜ x ∈ ∂ ˜Ω (120) ρ i (˜ x ,
0) = q i (˜ x ) for ˜ x ∈ ˜Ω , (121)where z i is the valence and D i is the diffusion coefficient for the ion specie i . The electricpotential φ (˜ x , t ) in ˜Ω is solution of the Neumann problem for the Poisson equation∆ φ (˜ x , t ) = − eε r ε ( ρ p (˜ x ) − ρ m (˜ x )) for ˜ x ∈ ˜Ω (122) ∂φ ( x , t ) ∂n = − ˜ σ (˜ x , t ) for ˜ x ∈ ∂ ˜Ω , where ˜ σ (˜ x , t ) is the surface charge density on the boundary ∂ ˜Ω. At steady-state, (119)gives ρ i (˜ x ) = ρ i, exp (cid:18) − z i eφ (˜ x ) k B T (cid:19) for i ∈ { p , m } , (123)where ρ i, is obtained from no-flux boundary condition (120), thus ρ i (˜ x ) = N i exp (cid:18) − z i eφ (˜ x ) k B T (cid:19)Z ˜Ω exp (cid:18) − z i eφ ( s ) k B T (cid:19) d s for i ∈ { p , m } . (124)27sing the non-dimensionalized potential ˜ u (˜ x ) = e φ (˜ x ) k B T , equation (123) becomes ρ i (˜ x ) = N i e − z i ˜ u (˜ x ) Z ˜Ω e − z i ˜ u ( s ) d s for i ∈ { p , m } . (125)Using (122) and (125), we obtain − ∆˜ u (˜ x ) = l B N p e − ˜ u (˜ x ) R ˜Ω e − ˜ u ( s ) d s − l B N m e ˜ u (˜ x ) R ˜Ω e ˜ u ( s ) d s in ˜Ω (126) ∂u (˜ x ) ∂n = − ( N p − N m ) | ∂ ˜Ω | l B on ∂ ˜Ω , where l B is the Bjerrum length. Using x = ˜ x R c and ˜ u (˜ x ) = u ( x ) where R c is the cuspcurvature radius, (126) becomes − ∆ u ( x ) = l B N p e − u ( x ) R c R Ω e − u ( s ) d s − l B N m eu ( x ) R c R Ω eu ( s ) d s in Ω (127) ∂u ( x ) ∂n = − l B ( N p − N m ) R c | ∂ Ω | on ∂ Ω . The small parameter is ζ = N m N p ≪ u ( x ) is u ( x ) = u ( x ) + ζ u ( x ) + · · · (128)Using (128) in (127), in small ζ limit, we have − ∆ u ( x ) = l B N p e − u ( x ) R c R Ω e − u ( s ) d s in Ω (129) ∂u ( x ) ∂n = − l B N p R c | ∂ Ω | on ∂ Ω , and in Ω ∆ u ( x ) = l B N p R c e − u ( x ) R Ω e − u ( s ) d s u ( x ) − R Ω u ( s ) e − u ( s ) d s R Ω e − u ( s ) d s ! + eu ( x ) R Ω eu ( s ) d s ! ∂u ( x ) ∂n = l B N p R c | ∂ Ω | on ∂ Ω , (130) ζ tends to zero (small charge limit) gives v ( x ), and thus we conclude that a smallamount of negative charges does not perturb the distribution of the total excess of positivecharge in the bulk. The numerical procedure
Numerical solutions were constructed by the Comsol Multiphysics 5.0 (BVP problems),Maple 2015 (Shooting problems) and Matlab R2015 (Conformal mapping). The boundaryvalue problems in 1D, 2D, and 3D were solved by the finite elements method in theComsol ’Mathematics’ package. We used an adaptive mesh refinement to ensure numericalconvergence for large value of the parameter λ .We solved the PDEs by the shooting procedure for boundary value problems usingRunge-Kutta 4 method, as well as solvers from Maple packages. References [1] Abramowitz, M., and I.A. Stegun,
Handbook of Mathematical Functions . Dover,(1972).[2] Bezanilla, F. ”How membrane proteins sense voltage.”
Nat Rev Mol Cell Biol. ,pp.323–332 (2008).[3] Bourne, J.N. and K.M. Harris, ”Balancing structure and function at hippocampaldendritic spines.” Annu Rev Neurosci. , , pp.47–67 (2008).[4] Cartailler, J., Z. Schuss, and D. Holcman, ”Analysis of the Poisson-Nernst-Planckequation in a ball for modeling the voltage-current relation in neurobiological mi-crodomains,” Physica D: Nonlinear Phenomena, , pp.39–48 (2016).[5] Courant, R. and D. Hilbert,
Methods of Mathematical Physics, v.2 , Wiley Inter-science, NY.(1989).[6] Debye, P. and E. H¨uckel, ”Zur Theorie der Elektrolyte. I. Gefrierpunktserniedrigungund verwandte Erscheinungen, Physikalische Zeitschrift,” , (9), pp.185–206 (1923).[7] Eisenberg, R.S., M.M. Klosek, and Z. Schuss, ”Diffusion as a chemical reaction:Stochastic trajectories between fixed concentrations.” J. Chem. Phys. , (4), 1767–1780 (1995).[8] Eisenberg, R.S., ”From structure to function in open ionic channels.” J. MembraneBiol. , , pp.1–24 (1998).[9] Eisenberg, R.S. ”Ionic channels in biological membranes. Electrostatic analysis of anatural nanotube. ” Contemp. Phys. , (6), pp.447–466 (1998).2910] Holcman, D. and Z. Schuss, ”Brownian motion in dire straits.” SIAM. J. on Multi-scale Modeling and Simulation (4), pp.1204–1231 (2012).[11] Holcman, D., N. Hoze, and Z. Schuss, ”Narrow escape through a funnel and effectivediffusion on a crowded membrane, Phys. Rev. E , , 021906 (2011); Erratum Phys.Rev. E , 039903 (2012).[12] Holcman, D. and Z. Schuss, ”Brownian motion in dire straits,” Multiscale Modeling& Simulation , , (4), pp.1204–1231 (2012).[13] Holcman, D. and Z. Schuss, Stochastic Narrow Escape in Molecular an Cellular Bi-ology, Analysis and Applications,
Springer Verlag, NY (2015).[14] Holcman, D. and R. Yuste, ”The new nanophysiology: regulation of ionic flow inneuronal subcompartments.”
Nature Reviews Neuroscience , pp.685–692 (2015).[15] Horn R., B. Roux, J. Aqvist, ”Permeation redux: thermodynamics and kinetics of ionmovement through potassium channels.” Biophys J. , (9), pp.185918–63 (2014).[16] Mamonov, A., R. Coalson, A. Nitzan, M. Kurnikova, ”The role of the dielectricbarrier in narrow biological channels: a novel composite approach to modeling singlechannel currents”, Biophys. J. , pp.3646–3661 (2003).[17] Morrison, J.A. and J.A. Lewis, ”Charge Singularity at the Corner of a Flat Plate,” SIAM Journal on Applied Mathematics , (2), pp.23-3-250 (1976).[18] Perry, D., D. Momotenko, R.A. Lazenby, M. Kang, and P.R. Unwin ”Characterizationof Nanopipettes”, Anal. Chem. , , (10), pp 5523–30, (2016).[19] Qian, N., T.J. Sejnowski, ”An electro-diffusion model for computing membrane po-tentials and ionic concentrations in branching dendrites, spines and axons.” Biol.Cybern. , pp.1–15 (1989).[20] Rall, W. ”Cable Theory for Dendritic Neurons. In Methods in Neuronal Modeling:from Synapses to Networks,” C. Koch, and I. Segev, eds. Cambridge, Mass., MITPres , pp.9–63 (1989).[21] Roux, B. and M. Karplus, ”Ion transport in the gramicidin channel: free energyof the solvated right-handed dimer in a model membrane,”
J. American ChemicalSociety , , pp.3250–3262 (1993).[22] Ruiz, F.J.G., A. Godoy, F. Gamiz, C. Sampedro, L. Donetti, ”A ComprehensiveStudy of the Corner Effects in Pi-Gate MOSFETs Including Quantum Effects,” IEEETransactions on Electron Devices, 54, (12), pp. 3369–3377 (2007).3023] Savtchenko, L.P., N. Kulahin, S.M. Korogod, D.A. Rusakov, ”Electric fields ofsynaptic currents could influence diffusion of charged neurotransmitter molecules.” Synapse , (4), pp.270–278 (2004).[24] Schuss, Z., B. Nadler, R.S. Eisenberg, ”Derivation of Poisson and Nernst-Planckequations in a bath and channel from a molecular model,” Physical Review E , SIAM of AppliedMathematics , ,3, pp.949–968 (2009).[26] Sparreboom, W., A. van den Berg and J.C.T. Eijkel, ”Principles and applications ofnanofluidic transport”, Nat. Nano. , , 713–720 (2009).[27] Yuste, R. Dendritic Spines . The MIT Press, Cambridge, MA (2010).[28] Lindsay, A., A. Bernoff, and M. Ward, ”First Passage Statistics for the Capture ofa Brownian Particle by a Structured Spherical Target with Multiple Surface Traps,”
SIAM J. Multiscale Modeling and Simulation, (2016).[29] Pillay, S., A. Peirce, T. Kolokolnikov, and M. Ward, ”An Asymptotic Analysis of theMean First Passage Time for Narrow Escape Problems: Part I: Two-DimensionalDomains,”
SIAM Multiscale Modeling and Simulation, (3), pp. 803-835 (2010).[30] Delgado, M. and D. Coombs, ”Conditional Mean First Passage Times to Small Trapsin a 3-D Domain with a Sticky Boundary, SIAM J. Multiscale Analysis and Simula-tion,13