Modeling of Ion and Water Transport in the Biological Nanopore ClyA

In recent years, the protein nanopore cytolysin A (ClyA) has become a valuable tool for the detection, characterization and quantification of biomarkers, proteins and nucleic acids at the single-molecule level. Despite this extensive experimental utilization, a comprehensive computational study of ion and water transport through ClyA is currently lacking. Such a study yields a wealth of information on the electrolytic conditions inside the pore and on the scale the electrophoretic forces that drive molecular transport. To this end we have built a computationally efficient continuum model of ClyA which, together with an extended version of Poison-Nernst-Planck-Navier-Stokes (ePNP-NS) equations, faithfully reproduces its ionic conductance over a wide range of salt concentrations. These ePNP-NS equations aim to tackle the shortcomings of the traditional PNP-NS models by self-consistently taking into account the influence of both the ionic strength and the nanoscopic scale of the pore on all relevant electrolyte properties. In this study, we give both a detailed description of our ePNP-NS model and apply it to the ClyA nanopore. This enabled us to gain a deeper insight into the influence of ionic strength and applied voltage on the ionic conductance through ClyA and a plethora of quantities difficult to assess experimentally. The latter includes the cation and anion concentrations inside the pore, the shape of the electrostatic potential landscape and the magnitude of the electro-osmotic flow. Our work shows that continuum models of biological nanopores—if the appropriate corrections are applied—can make both qualitatively and quantitatively meaningful predictions that could be valuable tool to aid in both the design and interpretation of nanopore experiments.


Introduction
The transport of ions and molecules through nanoscale geometries is a field of intense study that uses both experimental, theoretical and computational methods. [1][2][3][4][5][6] One of the primary driving forces behind this research is the development of nanopores as label-free, stochastic sensors at the ultimate analytical limit (i.e., single molecule). 7-10 Such detectors have applications ranging from the analysis of biopolymers such as DNA [11][12][13][14][15][16][17][18] or proteins, [19][20][21][22][23] to the detection and quantification of biomarkers, [24][25][26][27][28][29][30] to the fundamental study of chemical or enzymatic reactions at the single molecular level. 22,[31][32][33][34] Nanopores are typically operated in the resistive-pulse mode, where the fluctuations of their ionic conductance are monitored over time. [7][8][9]35 Experimentally, this is achieved by placing the nanopore between two electrolyte compartments and applying a constant DC (or AC) voltage across them. Due to the high resistance of the nanopore, virtually the full potential change occurs within (and around) the pore, resulting in a strong electric field (10 6 -10 7 mV · nm −1 ) that electrophoretically drives ions and water molecules through it. [36][37][38][39] Hence, analyte molecules such as DNA or proteins are driven towards, and often through, the nanopore by a combination of Coulombic (electrophoretic) and hydrodynamic (electro-osmotic) forces. 36,[40][41][42] If successful, a translocation event is observed as a temporal fluctuation in the ionic conductance of the pore that serves as a unique molecular 'fingerprint' with which the molecule can be identified. 43 Because the frequency, magnitude, duration and even noise levels of these events depend on the properties of both the analyte molecule and the nanopore itself, they are notoriously difficult to interpret unambiguously without a full understanding of the nanofluidic phenomena that underlie them.
However, the large computational cost of simulating a complete biological nanopore system (100K-1M atoms) for hundreds of nanoseconds still necessitates the use of supercomputers. 46,48 The PNP(-NS) equations, on the other hand, are of particular interest due to their low computational demands and analytical tractability. In a continuum approach, the simulated system is subdivided in several 'structureless' domains, the behavior of which is parameterized by material properties such as relative permittivity, diffusion coefficient, electrophoretic mobility, viscosity and density. Because these properties can only emerge from the collective behavior or interactions between small groups of atoms (i.e., the mean-field approximation), great care must be taken when using them to compute fluxes and fields at the nanoscale, where computational elements may only contain a few molecules. 65,66 Nevertheless, the PNP equations have been used extensively for the simulation of ion channels, 53,67,68 biological nanopores 58,63,69,70 and their solid-state counterparts 39,71-73 -often with excellent qualitative, if not quantitative results. 3,4,6 To remedy the shortcomings of PNP and NS theory, a number of modifications have been proposed over the years. These include, among others, 1) steric ion-ion interactions, 2) the effect of protein-ion/water interactions on their motility (i.e., diffusivity and electrophoretic mobility), 3) the concentration dependencies of ion motility, and solvent relative permittivity, viscosity and density.
The steric ion-ion interactions can be accounted for by computing the excess in chemical potential (µ ex i ) resulting from the finite size of the ions. 61 The interaction of ions or small molecules such as water with the heavy atoms of proteins or DNA results in a strong reduction of their motility, as observed in MD simulations. 79,80 Since these effects happen only at distances ≤1 nm, they can usually be neglected for macroscopic simulations. However, in small nanopores (≤ 10 nm radius), they comprise a significant fraction of the total nanopore radius and hence must be taken into account. 54,58,63,81 In continuum simulations, this can be achieved with the use pf positional-dependent ion diffusion coefficients. An example implementation is the 'soft-repulsion PNP' developed by Simakov and Kurnikova,63,70 who used it to predict the ionic conductance of the αhemolysin nanopore. Similar reductions in ion diffusion coefficients have been proposed to improve PNP theory's estimations of the ionic conductance of ion channels. 67,68 The motility of water molecules is expressed by the NS equations as the fluid's viscosity. Hence, as also observed in MD simulations for water molecules near proteins 80 and confined in hydrophilic nanopores, 82-84 the water-solid interaction leads to a viscosity several times higher compared to the bulk values. Note that this is valid for hydrophilic interfaces only, as the lack of interaction with hydrophobic interfaces, such as carbon nanotubes, leads to a lower viscosity. 85 It is well known that the self-diffusion coefficient D i and electrophoretic mobility µ i of an ion i depends on the local concentrations of all the ions in the electrolyte. 86 Their values typically decrease with increasing salt concentration, and should not be treated as constants.
Moreover, even though the well-known Nernst-Einstein (NE) relation µ i = D i /k B T is strictly speaking only valid at infinite dilution and a good approximation at low concentrations (<10 mM), it significantly overestimates the ionic mobility at higher salt concentrations. [86][87][88][89] In an empirical approach, Baldessari and Santiago formulated an ionic-strength dependency of the ionic mobility based on the activity coefficient of the salt 60 and showed excellent correspondence between the experimental and simulated ionic conductance of long nanochannels over a wide concentration range. 90 Alternatively, Burger et al. used a microscopic lattice-based model to derive a set of PNP equations with non-linear, ion density-dependent mobilities and diffusion coefficients that provided significantly more realistic results for ion channels. 91 Note that other electrolyte properties, such as its viscosity, 92 density 92 and relative permittivity, 93 also significantly affect the ion and water flux. To better compute the charge flux in ion channels, Chen derived a new PNP framework 94 that includes water-ion interactions in the form of a concentration-dependent relative permittivity and an additional ion-water interaction energy term.
To the best of our knowledge, no attempt has been made to consolidate all of the corrections discussed above into a single framework. Hence, we propose an extended set of PNP-NS (ePNP-NS) equations, which improves the predictive power of the PNP-NS equations at the nanoscale and beyond infinite dilution. Our ePNP-NS framework takes into account the finite size of the ions using a size-modified PNP theory, 78

Mathematical model
The use of continuum or mean-field representations for both the nanopore and the electrolyte enables us to efficiently compute the steady-state ion and water fluxes under almost any condition. The dynamic behavior of our complete system is described by the coupled Poisson, Nernst-Planck and Navier-Stokes (PNP-NS) equations, a well-known set of partial differential equations that describe the electrostatic field, the total ionic flux and the fluid flow, respectively. 61,64,71 In this section we will describe all the components of the simulation, i.e., the 2D-axisymmetric nanopore geometry and the system of governing equations. The 2D-axisymmetric geometry was derived directly from the all-atom model by computing the average inner and outer radii along the longitudinal axis of the pore, and hence closely follows the outline of a 30°wedge out of the homology model. (d) The fixed space charge density (ρ f pore ) map of ClyA-AS, obtained by Gaussian projection of each atom's partial charge onto a 2D plane (see methods for details). (e+f) The 2D-axisymmetric simulation geometry of ClyA (grey) embedded in a lipid bilayer (green) and surrounded by a spherical water reservoir (blue). Note that all electrolyte parameters depend on the local average ion concentration c = 1 n n i c i and that some are also influenced by the distance from the nanopore wall d.

Model geometry
2D-axisymmetric model of ClyA ClyA is a relatively large protein nanopore that selfassembles on lipid bilayers to form 14 nm long hydrophilic channels. The interior of the pore can be divided into roughly two cylindrical compartments (Fig. 1a): the cis lumen (≈6 nm diameter, ≈10 nm height), and the trans constriction (≈3.3 nm diameter, ≈4 nm height).
Because ClyA consists of 12 identical subunits (Fig. 1b), it exhibits a high degree of radial symmetry, a geometrical feature that can be exploited to obtain meaningful results at a much lower computational cost. 58,64,71 However, this requires the reduction of the full 3D atomic structure and charge distribution to a realistic 2D-axisymmetric model. To this end, we con- that are impermeable to ions and water.

Governing equations
In an attempt to improve upon the quantitative accuracy of the PNP-NS equations for nanopore simulations, we developed an extended version of these equations (ePNP-NS) and implemented it in the commercial finite element solver COMSOL Multiphysics (v5.4, www.comsol.com). Our ePNP-NS equations self-consistently take into account 1) the finite size of the ions, 76,78 2) the reduction of ion and water motility close to the nanopore walls, 54,58,79,80,83 and 3) the concentration dependency of ion diffusion coefficients and electrophoretic mobilities, as well as electrolyte viscosity, density and relative permittivity. 87,92,93 Most of these corrections make use of using empirical functions that were fitted to experimental data (Tabs. 1 and S2).
The pore's fixed charge distribution, ρ f pore , (see Eq. 18) was derived directly from the full atom model of ClyA-AS. The ionic charge density in the fluid is given by with F Faraday's constant (96 485.33 C · mol −1 ), and c i the ion concentration and z i ion charge number of ion i.
To account for the concentration dependence of the electrolyte's relative permittivity, we used the expression with c = 1 n n i c i the average ion concentration, ε 0 r the relative permittivity at infinite dilution and ε c r ( c ) a concentration dependent empirical function parameterized with experimental data.
Ionic flux. The total ionic flux J i of each ion i is given by the size-modified Nernst-Planck equation, 78 and can be expressed as the sum of diffusive, electrophoretic, convective and steric fluxes where and at steady state The reduction of the ionic motility at increasing salt concentrations and in proximity Ion self-diffusion coefficient Ion electrophoretic mobility to the nanopore walls was implemented self-consistently by replacing D i and µ i with the where Based on the observation that the diffusivity of nanometer-to micrometer-sized particles reduces significantly when confined in pores and slits of comparable dimensions, 42,[110][111][112][113] Simakov et al. 63 and Pederson et al. 58 reduced the ion motilities inside the pore as a function of the ratio between the ion and the nanopore radii. We chose not to include this correction into our model, as extrapolating its applicability for ions with a hydrodynamic radii comparable to size of the solvent molecules is questionable. 111,114 Fluid flow. The fluid flow and pressure are given by the Navier-Stokes equations for incompressible fluids 115 where together with the continuity equations for the fluid density u · ∇ = 0 , (11) and the fluid velocity with u the fluid velocity, the fluid density, σ ij the hydrodynamic stress tensor, η the viscosity and p the pressure. The external body force density F that acts on the fluid is given by with E = −∇ϕ the electric field vector.
As with the previous equations, we introduced a concentration dependency and wall distance dependencies for η and a concentration dependency for the by replacing their constant values by where η 0 and 0 are the values at infinite dilution (i.e., pure water). The empirical functions η c ( c ), c ( c ) and η w (d) were parameterized via fitting to experimental 92 and molecular dynamics 80 data obtained from literature.

Boundary conditions
The reservoir boundaries were set up, with Dirichlet conditions, to act as electrodes: the cis side was grounded (ϕ = 0) and a fixed, but changeable, bias potential was applied

Note on the concentration dependencies.
All concentration dependent parameters use the local ionic strength rather than their individual ion concentrations. Though valid for electroneutral bulk solutions, this approximation no longer holds inside the electrical double layer (i.e., near charged surfaces or inside small nanopores), where local electroneutrality is violated. The main reasons for nevertheless making this simplification are the lack of non-bulk experimental data and the absence of a tractable analytical model. Furthermore, we will see that the current implementation of our concentration dependent functions will lead to an excellent agreement with the experimental data in all but the most extreme cases, justifying our choice a posteriori.

PNP-NS vs. ePNP-NS
The ePNP-NS equations revert into the regular PNP-NS equations disabling the steric flux (β = 0) and by setting all concentration and wall distance functions to unity (ε r = ε 0 r ,

Results and Discussion
The current-voltage (IV) relationships of many nanopores, ClyA included, often deviate significantly from Ohm's law. This is because the ionic flux arises from a complex interplay between the pore's geometry (e.g., size, shape and charge distribution), the properties of the surrounding electrolyte (e.g., salt concentration viscosity and relative permittivity) and the externally applied conditions (e.g., bias voltage, temperature and pressure). The ability of a computational model to quantitatively predict the ionic current of a nanopore over a wide range of bias voltages and salt concentrations strongly indicates that it captures the essential physics governing the nanofluidic transport. Hence, to validate our model, we experimentally measured the single channel ionic conductance of ClyA at a wide range of experimentally relevant salt concentrations and bias voltages. We compared these experimental data with the simulated ionic transport properties in terms of current, conductance, rectification and ion selectivity, of both the classical PNP-NS and the newly developed ePNP-NS equations.
After validation, we proceed with describing the influence of the bias voltage (V b ) and bulk salt concentration (c s ) on the local ion and charge distribution and inside the pore. This is followed by a characterization of the electrostatic potential and the electrostatic energy landscape within ClyA for both cations and anions. We will conclude this section by discussing the properties of the electro-osmotic flow.

Transport of Ions Through ClyA
Ionic current and conductance. The ability of our model to reproduce the ionic current of a biological nanopore over a wide range of experimentally relevant conditions (between V b = −150 to +150 mV and for c s = 0.05, 0.15, 0.5, 1 and 3 M NaCl) can be seen in Fig. 2a.
Here, we compare IV relationships of ClyA-AS as measured experimentally ('expt.'), simulated using our 2D-axisymmetric model ('PNP-NS' and 'ePNP-NS') and naively analytically estimated ('bulk') using a resistor model of the pore 96,116 (Eqs. S18 and S19). Whereas the   The ability of a nanopore to conduct ions can be best expressed by its conductance: The cross-sections of the ionic conductance as a function of concentration at high positive and negative bias voltages (Fig. 2c), serve to demonstrate the differences between these respective regimes. The slope of the log-log plot at +150 mV is highly linear, suggesting that only a single mechanism of ionic conduction is at work. As touched upon above, both the bulk model and the ePNP-NS equations manage to capture the conductance at high ionic strengths, while only the latter performs well at low ionic strengths. The predictions made with the PNP-NS equations overestimate the conductance over the entire concentration range, but they do converge at very low ionic strenghts. At V b = −150 mV, the log-log plot consists of two linear segments with a transition zone at c s ≈ 0.15 M. This behavior is captured qualitatively by both simulation methods, but a perfect quantitative match is only found for the ePNP-NS equations. The bulk model exhibits the same single-sloped trend as The difference in ionic conduction at opposing bias voltages is also known as ionic current ICR is a phenomenon often observed in nanopores that are both charged, and contain a degree of geometrical asymmetry along the central axis of the pore. 5,72,117 As can be seen in Fig. S3, ClyA exhibits a strong degree of rectification, which is to be expected given its high negative charge (−72 e at pH 7.5) and different cis (≈3.3 nm) and trans (≈6 nm) entry diameters (Fig. 1a). We found α to increase monotoneously with the bias voltage magnitude, at least over the investigated range. We found the dependence of α on the ionic strength not to be monotonous, but rather rising rapidly to a peak value at ≈0.15 M, followed by a gradual decline towards unity at saturating salt concentrations.
The results and comparisons discussed above indicate that ClyA's conductivity is dominated by the bulk electrolyte conductivity above physiological salt concentrations (c s > 0.15 M). The breakdown of this simple dependency at lower ionic strengths is particularly evident at negative bias voltages and is likely caused by the overlapping of the electrical double layer inside the pore. This effectively excludes the co-ions-Clin the case-from the interior of the pore, preventing them from contributing to total ionic conductance and resulting in a conductance dominated by the surface charge. 118 The presence of only a single ion type inside the pore may also offer an explanation as to why the ePNP-NS equations falter at low ionic strengths. Because our ionic mobilities are derived from bulk ionic conductances, i.e., for unconfined ions in a locally electroneutral environment, it is likely that our mobility model begins to break down under these conditions. 119 Another cause of the discrepancies could be a change in the shape or diameter of the nanopore at low salt concentrations, which cannot be captured by our simulation due to the static nature of its geometry and charge distributions. Nevertheless, our simplified 2D-axisymmetric model, in conjunction with the ePNP-NS equations, is able to accurately predict the ionic current flowing through ClyA for a wide range of experimentally relevant ionic strengths and bias voltages. This suggests that we managed to capture the essential physical phenomena that drive the ion and water transport through the nanopore both qualitatively and quantitatively. Hence, we expect the distribution of the resulting properties (e.g., ion concentrations, ion fluxes, electric field and water velocity) to closely correspond to their real values.
Cation selectivity. The ion selectivity of a nanopore determines the preference with which it transports one ion type over the other. Experimentally, it is often determined by placing the pore in a salt gradient (i.e., different salt concentrations in the cis and trans reservoirs) and measuring the potential at which the nanopore current is zero (reversal potential, V r ). 96,101 The Goldman-Hodkin-Katz (GHK) equation can then be used to convert V r into the permeability ratio P Na + = G Na + /G Cl − . Here, we represent the ClyA's ion selectivity (Figs. 2d and 2e) by the fraction of the total current that is carried by Na + ions: the apparent Na + transport number t Na + = G Na + /(G Na + + G Cl − ) = P Na + /(P Na + + 1).
Due its negatively charged interior, we found ClyA to be cation selective (i.e., t Na + > 0. 5) for all investigated voltages up to a bulk salt concentration of c s ≈ 2 M NaCl (0.5 contour line in Fig. 2d). Above this concentraion, t Na + falls to a minimum of value of 0.45 at c s ≈ 5 M, which is still ≈1.27 times its value in the bulk electrolyte (0.35). This indicates that-even at saturation-ClyA remains preferential towards cations. Below 2 M, the ion selectivity increases logarithmically with decreasing salt concentrations, but it also becomes more sensitive to the direction and magnitude of the electric field: with negative bias voltages yield higher ion selectivities (Fig. 2e). For example, to reach a selectivity of t Na + ≈ 0.9, the salt concentration must fall to 0.05 M at +150 mV and to only 0.125 M at −150 mV.
Using the reversal potential method, Franceschini et al. 101 found ClyA's ion selectivity to be t Na + = 0.66 (P Na + = 1.9). In our ePNP-NS simultation, this corresponds the selectivity at c s = 0.5 M (at V b = 0 mV), a value that lies in between the cis (1 M, t Na + = 0.57, concentrations. These two effects should not be ignored as they contribute significantly to the nanopore's total conductance. Furthermore, because the ion selectivity depends strongly on the ionic strength, the measured reversal potential will necessarily be influenced by the chosen salt gradient and represent the selectivity at an undetermined intermediate concentration.  We also observed a stark difference between their sensitivities to the applied bias voltage, particularly at low salt concentrations (Fig. 3a,  of the pore lumen (d ≥ 0.5 nm), but with pockets of highly charged and alternating positive and negative charge densities close to the nanopore wall (Fig. 3e, rightmost panel). This sharp confinement is shown clearly by the radial density profiles (Fig. 3f) drawn through the constriction (z = −0.3 nm, purple triangles) and the lumen (z = 5 nm, green triangles).  , the negative electrostatic potential extends significantly inside the lumen (1.6 < z < 12.25 nm), and even more so inside the trans constriction (1.85 < z < 1.6 nm). For the former, localized influential negative 'hotspots' can be found in the middle (4 < z < 6 nm) and at the cis entry (10 < z < 12 nm). (c) Radial average of the equilibrium electrostatic potential along the length of the pore ( ϕ rad ) for the same concentrations as in (b). Even though the lumen of the pore is almost fully screened for cs > 0.5 M, the constriction still retains some of its negative influence even at 5 M.

Electrostatic Potential and Energy
The electrostatic potential, or rather the change thereof, is one of the primary driving forces within a nanopore. Typically, the potential can be split into an intrinsic, 'equilibrium' component and an external, 'non-equilibrium' contribution. The latter is also known as the bias voltage applied between the trans and the cis reservoirs, and most of it changes across the short distance of the nanopore itself. For a 14 nm-long nanopore such as ClyA, a potential drop of 10 to 100 mV can easily create an electric field of 10 6 -10 7 V · m −1 inside the pore.
The intrinsic electrostatic potential results from the charged residues embedded in the interior walls of the protein and can significantly influence the transport of ions and water molecules through the pore. 46,48,57,120 In the following section we aim to describe the most salient features of this potential and the extent of its influence over the entire investigated concentration range (Fig. 4). Moreover, because electromigration is the primary contributor to the ionic current, a proper understanding of the electrostatic landscape-and the energy barriers faced by the ions traversing the pore-should provide a more quantitative explanation for the origin of ClyA's rectification, ion selectivity and asymmetric ion concentrations.
Distribution of the equilibrium electrostatic potential. The electrostatic potential at equilibrium (ϕ 0 , i.e., at V b = 0 mV) reveals the effect of ClyA's fixed charges on the potential inside the pore (Fig. 4b). Due to electric screening by the mobile charge carriers in the electrolyte, however, the extent of their influence strongly depends on the bulk ionic strength. The contour plot cross-sections of ϕ 0 for c s = 0.005, 0.05, 0.15, 0.5 and 5 M (Fig. 4b) and their corresponding radial averages (Fig. 4c) demonstrate this effect aptly. The radial average ( ϕ 0 rad ) represents the mean value along the longitudunal axis of the pore and can be computed using where R(z) is taken as the pore radius for −1.85 ≤ z ≤ 12.25 nm, 2 nm for z < −1.85 nm and 4 nm for z > 12.25 nm. Starting from the cis entry (z ≈ 10 nm), the electrostatic potential is dominated by the acidic residues D114, D121 and D122, resulting in a rapid reduction of ϕ 0 rad upon entering the pore. Next, ϕ 0 rad slowly descreases up until the middle of the lumen (z ≈ 5 nm), where the next set of negative residues, namely E53, E57 and D64, lower it even further. After a brief increase, ϕ 0 rad attains its maximum amplitude inside the trans constriction (z ≈ 0 nm) due to the close proximity of the amino acids E7, E11, E18, D21 and D25, and then quickly falls to 0 inside the trans reservoir. Non-equilibrium electrostatic energy at +150 and −150 mV. To link back the observed ionic conductance properties to the electrostatic potential, we computed the radially averaged electrostatic energy for a monovalent ion, U E,i rad = z i e ϕ rad , at V b = +150 −150mV for the entire range of simulated ionic strengths (Fig. 5a). The resulting energy plot represents the energy landscape-filled with barriers (hills) or traps (valleys)-that a positive or negative ion must traverse in order to contribute positively to (i.e., increase) the ionic current.  At positive bias voltages, cations traverse the pore from trans to cis (Fig. 5a, first plot).
Upon entering the negatively charged constriction, their electrostatic energy drops dramatically, followed by a relatively flat section with a small barrier for entry in the lumen at At negative voltages, cations move through the pore from cis to trans, with a slow and continuous drop of the electrostatic energy throughout the lumen of the pore up until the constriction (Fig. 5a, third plot). This results in the efficient removal of cations from the pore lumen, and explains the lower Na + concentration observed at positive voltages (Fig. 3a).
To fully exit from the pore, however, cations must overcome a large energy barrier, which reduces the nanopore's ability to conduct cations compared to positive potentials and hence contributes to the ion current rectification. The situation for anions at negative bias voltages (i.e., travelling from trans to cis) is very different (Fig. 5a, fourth plot). Their ability to even enter the pore is severely hampered by an energy barrier of a few k B T at the trans constriction. Any anions that do cross this barrier, and those still present in the lumen of ClyA, will be rapidly move towards the cis entry and exit from the pore due to a continuous drop of electrostatic energy. This effectively depletes the entire lumen of anions, which can be observed from the much lower Clconcentrations at negative voltages (see Fig. 3a).
Concentration and voltage dependencies of the energy barrier at the constriction. Many biological nanopores contain constrictions that play crucial roles in shaping their ionic conductance properties. 14,28,101 The reason for this is two-fold, 1) the narrowest part dominates the overall resistance of the pore and 2) confinement of charged residues results in much larger electrostatic energy barriers. With its highly negatively charged trans constriction, ClyA's affinity for transport of anions is diminished and that for cations is enhanced compared to bulk, even at high ionic strengths (Fig. 2e). 96 To further elucidate the significance of the trans electrostatic barrier (∆E B,i ), we quantified its height at positive and negative voltages as a function of the salt concentration (Fig. 5b).
Due to the tilting of the nergy landscape, the application of a bias voltage lowers the

Transport of Water Through ClyA
The charged nature of many nanopores gives rise to a net flux of water through the pore, called the electro-osmotic flow (EOF). 37,82,121 The EOF not only contributes significantly to the ionic current, but the magnitude of the viscous drag force it exerts is often of the same order as the Coulombic electrophoretic force (EPF). Hence, it strongly influences the capture and translocation of biomolecules including nucleic acids, 36,47,122 peptides, 28,123,124 and proteins. 25,27,30,[96][97][98][99]125 Because the drag exerted by the EOF depends primarily on the size and shape of the biomolecule of interest and not on its charge, 125 it can be harnassed to capture molecules even against the electric field. 25 The EOF is a consequence of interaction between the fixed charges on the nanopore walls and mobile charges in the electrolyte. It can be described by two closely related mechanisms: 1) the excess transport of the hydration shell water molecules in one direction due to the pore's ion selectivity, and 2) the viscous drag exerted by the unidirectional movement of the electrical double layer inside the pore.
The first mechanism likely dominates in pores with a diameter close to that of the hydrated ions (≤1 nm) such as αHL or FraC, 28,124 while the second is expected to be stronger for larger pores (>1 nm), such as ClyA 25,125 or most solid-state nanopores. 37,39 In our simulation, the EOF follows the coupling of the Navier-Stokes and the Poisson-Nernst-Planck equations through a volume force term (Eq. 13), which dictates that the electric field exerts a net force on the fluid if it contains a net charge density-as is the case for the electrical double layer lining the walls of ClyA (Fig. 3d). In the following section we aim to quantitatively and a qualitatively describe the properties of the water velocity inside ClyA and how it is influenced by the bulk ionic strength and the applied bias voltage (Fig. 6).
Direction, magnitude and distribution of the water velocity. As expected, given ClyA's negatively charged interior surface and the resulting positively charged electrical double layer, the direction of the net water flow inside ClyA follows the electric field, i.e., from the cis to the trans at negative bias voltages (Fig. 6a). This corresponds to the observations and analysis of single-molecule protein capture 96 and trapping 97,98,125 experiments using the ClyA-AS nanopore. Along the longitudinal axis (z) of the pore, the water velocity is governed by the conservation of mass, meaning it is lowest in the wide cis lumen and highest in the narrow trans constriction (Fig. 6a). For example, at V b = −100 mV and c s = 0.5 M the velocity at the center of the pore is ≈0.07 m · s −1 in the lumen and ≈0.21 m · s −1 in the constriction.
Along the radial axis (r), u has a parabolic profile with the highest value at the center of the pore and the lowest at the wall due to the no-slip boundary condition (Figs. 6b and 6c). This is the result of a self-induced pressure gradient caused by the expansion of the EOF as it exits the pore. 126 Influence of bulk ionic strength and bias voltage on the electro-osmotic conductance. In analogy to the ionic conductance, the amount of water transported by ClyA can be expressed by the electro-osmotic conductance G eo = Q eo /V b (Fig. 6d). Here, Q eo is the net volumetric flow rate of water through the pore and computed by integrating the water velocity across the reservoir boundary (Eq. 21). The strength of the EOF depends strongly and non-monotonically on the bulk ionic strength: G eo rapidly increases with ionic strength until a peak value is reached at c s ≈ 0.5 M, followed by a gradual logarithmic decline (Fig. 6d).
The sensitivity of the EOF to the magnitude and sign of the bias voltage is given by the electro-osmotic conductance rectification α eo (V ) = G eo (+V )/G eo (−V ) (Fig. 6e). For all voltages magnitudes, α eo shows a maximum at c s ≈ 0.045 M, after which it falls rapidly to reach unity (α eo = 1) at approximately c s ≈ 0.45 M. A minimum is then reached at ≈1 M, followed by a gradual approach towards unity at c s = 5 M.

Conclusions
We We have verified our approach by matching experimental results to a very high degree of accuracy and the model parameters where gauged on other experiments, ultimately leaving no degrees of freedom. This shows that a continuum approach to modeling biological nanopores is not only feasible but to a very high degree predictive. We are confident that such a model as presented in this paper, constitutes a highly practical tool that can help with 1) elucidating the link between ionic current observed during a nanopore experiment and the actual physical phenomenon, 2) describing the electrophoretic and electro-osmotic properties of any biological nanopore and 3) guiding the design of new variants of existing nanopores.
In the analysis of the operation of ClyA with our model we found that the ionic currents depend strongly on the details of the electric field structure within the pore which define potential barriers for the ion species. These barriers in turn are influenced by both the bulk ionic strength and the applied bias voltage. Nevertheless, despite its simplified 2Daxisymmetric nature, our model is able to obtain meaningful and predictive results in terms of the ionic current, the ion concentrations, the electrostatic potential and electro-osmotic flow velocity. Some inaccuracies still arise at low salt concentrations, but these can likely be resolved with an improved mobility model.
Although not explicitly investigated in the current work, there are no obvious reasons why the model framework is not directly applicable to other biological nanopores, solid-state nanopores and other nanoscale ionic transport problems. These systems can be treated analogously and their properties could be mapped out systematically using our ePNP-NS model, greatly enhancing the understanding and interpretation of experimental results.

Materials and Methods
ClyA-AS homology model. membrane builder 130 and equilibrated with NAMD 131 , as described in detail in ref. [ 132]. Axially symmetric geometry. The 2D-axisymmetric geometry of the ClyA-AS nanopore ( Fig. 1c) was derived directly from its full atom model by radially averaging the molecular density. Briefly, 50 sets of atomic coordinates were extracted from the final 5 ns of the coordinates of the 30 ns MD production run (i.e., every 100 fs) and aligned by minimizing the RMSD between their backbone atoms (VMD RMSD tool). Next, we computed and averaged the 3D-dimensional molecular density maps of all 50 structures on a 0.5 A resolution grid using the Gaussian function 104 where for each atom i, R i is its Van der Waals radius, is the distance of grid coordinates (x, y, z) from the atom center (x i , y i , z i ) and σ = 0.93 is a width factor. The resulting 3D density map was then radially averaged along the z-axis, relative to the center of the pore to obtain a 2D-axisymmetric density map. The contourline at 25 % density was used as the nanopore simulation geometry, after manual removal of overlapping and superfluous vertices to improve the quality of the final computational mesh.
Axially symmetric charge density. The 2D-axially symmetric charge distribution (Fig. 1d was also derived directly from the 50 sets of aligned nanopore coordinates) that were used for the geometry. Inspired by how charges are represented in the particle mesh Ewald (PME) method, 46 we computed the fixed charge distribution of the nanopore ρ f pore (r, z) by assuming that an atom i of partial charge δ i , at the location (x i , y i , z i ) in the full 3D atomistic pore model, contributes an amount δ i /2πr i to the partial charge at a point (r i , z i ) with r i = x 2 i + y 2 i in the averaged 2D-axisymmetric model. This effectively spreads the charge over all angles to achieve axial symmetry. We assumed a Gaussian distribution of the space charge density of each atom i around its respective 2D-axisymmetric coordinates (r i , z i ) is where R i is the atom radius, σ = 0.5 is the sharpness factor and e is the elementary charge.
To embed ρ f pore with sufficient detail, yet efficiently, into a numeric solver, the spatial coordinates were discretized with a grid spacing of 0.005 nm in the domain of ρ f pore and precomputed values were used during the solver runtime. All partial charges (at pH 7.5) and radii were taken from the CHARMM36 forcefield 133 and assigned using PROPKA 134 and PBD2PQR. 135 Computing electrophoretic mobilities. To obtain the concentration-dependent ionic mobility µ c i from fitted functions, it must first be derived from the salt's molar conductivity Λ and the ion's transport number t i before it can be fitted 86 where λ i (c) is the specific molar conductivity of ion i.
Computing the simulated ionic current and electro-osmotic flow rate. The simulated ionic current I sim at steady-state was computed by with z i the charge number and J i the total flux of each ion i across cis reservoir boundary S, F the Faraday constant (96 485 C · mol −1 ) andn the unit vector normal to S. Similarly, the volumetric flow rate, i.e., the volume of water passing through the pore per unit time, is given by ClyA expression and purification. ClyA-AS monomers were expressed, purified and oligomerized using methods described in detail elsewhere. 25 Recording of single-channel current-voltage curves. Experimental current-voltage curves where measured using single-channel electrophysiology, as detailed elsewhere. 25,35,96 Briefly, a black lipid bilayer was inside a ≈100 µm diameter aperture in a thin teflon film separating two buffered electrolyte compartments by painting with 1,2-diphytanoyl-snglycero-3-

Supporting Information Available
The supplementary info contains the Extended materials and methods, with details on the fitting of the electrolyte properties and the calculation of the pore averaged values, and the weak forms of the ePNP-NS equations. It also contains additional results, including a figures about ionic current rectification and the electro-osmotic pressure distribution inside the pore, and a table detailing the peak values of the radial electrostatic potential inside ClyA.