Transport of neutral and charged nanorods across varying-section channels

We study the dynamics of neutral and charged rods embedded in varying-section channels. By means of systematic approximations, we derive the dependence of the local diffusion coefficient on both the geometry and charge of the rods. This microscopic insight allows us to provide predictions for the permeability of varying-section channels to rods with diverse lengths, aspect ratios and charge. Our analysis shows that the dynamics of charged rods is sensitive to the geometry of the channel and that their transport can be controlled by tuning both the shape of the confining walls and the charge of the rod. Interestingly, we find that the channel permeability does not depend monotonically on the charge of the rod. This opens the possibility of a novel mechanism to separate charged rods.


I. INTRODUCTION
The transport of molecules, proteins and small particles across pores, channels and, in general, porous materials is of paramount relevance in several biological, environmental and technological scenarios [1]. In biology and technological applications, ions and molecules are transported across membrane pores [2][3][4], RNA is transported across the nuclear membrane [5][6][7] and the crowded environment of cell cytoplasm provides dynamical obstacles to the motion of organelles and vesicles [8]. Similarly, in environmental sciences, pollutants and plant nutrients spread across the porous matrix of rocks and soil. In this perspective, substantial effort has been dedicated to the synthesis [9,10] and characterization [11] of anisotropic nano-colloids and stiff filaments [12,13]. Finally, in recent years the development of experimental techniques such as DNA [14,15] and protein [16] sequencing, resistive-pulse sensing techniques [17][18][19][20] and chromatography [21][22][23] have been developed by exploiting transport across pores.
On the theoretical side, tackling these problems is typically complicated since it requires detailed numerical simulations that take into account the interactions between the suspended object and the confining walls. Recently, several groups have contributed developing the so-called Fick-Jacobs approximation [24][25][26][27][28][29]. This scheme relies on the assumption that when the net longitudinal velocity is slow enough, the transported object will explore the directions transverse to the motion with a probability distribution that is very close to the equilibrium one. Thanks to this approximation, the geometry of the pore as well as the interaction with the walls can be accounted for by the local equilibrium free energy. This approach has been shown to be applicable to the transport of colloids [30][31][32][33][34], ions [35][36][37][38][39], and polymers [40][41][42] just to mention a few among others (see recent review articles [29,43] for a more comprehensive list).
In this article we study the dynamics of charged rods confined in varying section channels. In particular, we propose an extension of the Fick-Jacobs approximation that allows us to account for both the geometry of the rod as well as for its net charge. From this perspective, our work goes beyond a recent contribution (see Ref. [44]) in which an approximated equation governing the transport of uncharged rods across varying section channels has been proposed. First, we test our approach against known experimental results [45] for neutrally charged rods. Our model retrieves both the local dependence of the diffusion coefficient as well as the dependence of the Mean First Passage Time (MFPT) profile. Second, we use our model to predict the dependence of the diffusion coefficient on the geometry of both the channel and the rod as well as on the charge of the rod. For neutrally charged rods, in agreement with Ref. [45], we find that the local diffusion coefficient maximizes at the channel bottleneck, i.e. where the rod is mainly oriented along the channel axes. Interestingly, the modulation of the local diffusion coefficient is maximized when the length of the rod major axis is comparable to the channel average section and it becomes vanishing small when these two length scales differ. Finally, our data show that such a dependence can be even amplified by the electrostatic interaction between the rod and the channel walls.

II. MODEL
We characterize the transport of rigid rods across a 2D varying-section channel whose half-section is given by where h 0 is the average section, h 1 is the amplitude of the section modulation and L 0 is the period of the channel. Accordingly, the 2D-Smoluchowski equation reads with where ρ is the probability distribution function, β −1 = k B T is the inverse thermal energy, k B is the Boltzmann constant, T is the absolute temperature. Moreover, D(θ) is the translational diffusion matrix that accounts for the fact that the diffusion along the major axis of the rod differs from that perpendicular to it, D θ is the rotational diffusion coefficient and W accounts for the both, the conservative forces acting on the rod and the confinement. In the following, we focus on the case in which the local radius of curvature is much larger than the length of the rod. In such a regime the channel walls can be approximated as locally parallel to the longitudinal axis. Accordingly we have that the potential W reads where φ(x, z, θ) accounts for the equilibrium conservative forces and f for the longitudinal force responsible for the transport along the channel. In Eqs. (7), (8) we assume that the aspect ratio between the long axis, 2L, and the minor axis, 2l, of the rod is such that L ≫ l and hence we can disregard the geometric corrections to Eqs. (7), (8) due to the finiteness of l.
We further assume that translation along the channel axis is slow enough such that ρ retains its equilibrium profiles along z and θ, which implies the absence of fluxes in z and θ, i.e.
Accordingly, we perform the standard Fick-Jacobs approximation [24][25][26][27][28] ρ(x, z, where is the local equilibrium free energy [24,25]. We remark that p(x) (see Eq. (11)) is proportional to the probability of finding the center of mass at a position x. In the frame of reference of the rod we havē where D 0 l L is the diffusion coefficient along the minor axis of the rod (of size l), e ⊥ , and D 0 is the diffusion coefficient along the major axis of the rod (of size L), e , and the off-diagonal terms are zero due to the axial symmetry of the rod [46]. Integrating Eq. (3) along z and θ and using Eqs.
Here, is the effective local diffusion coefficient. According to Fig. 1 we have e x = cos θe − sin θe ⊥ which leads to the following expression for the local diffusion coefficient, which, once substituted into Eq. (15), leads to (17) We remark that for l = L Eq. (17) reduces to D(x) = D 0 , which is the diffusion coefficient of spherical particles of size l [24,25]. Finally, integrating Eq. (2) along z and θ and using Eqs. (14) leads to (18) This is the extension of the Fick-Jacobs approximation to the case of rigid rods embedded in varying-section channels. The steady-state solution of Eq. (18) reads where J and Π are integration constants. For periodic boundary conditions and imposing the normalization of the probability we obtain Finally we define the dimensionless permeability as where f is the magnitude of the external force and J 0 = βD 0 f /L is the flux of a spherical colloid of radius l in a flat channel of half-width h 0 .

A. Neutral rods
At first we focus on neutral rods for which the potential W reads (23) and the local free energy reads where Θ(θ) is the Heaviside step function. Accordingly, the effective diffusion coefficient reads where L is the half-length of the major axis of the rod (see Fig. 1).

Diffusion
At first we compare our analytical predictions for the diffusion coefficient against experimental results. As already shown in Ref. [45], the diffusion coefficient is maximum at the channel bottleneck where the rod is mainly parallel to the axis of the channel. Fig. 2a shows that predictions of Eq. (26) qualitatively agree with the experimental data. In particular, the model properly captures the enhancement of the sensitivity of the local diffusion coefficient upon increasing the length of the rod. Indeed, we remark that for longer rods, for which the dependence of the diffusion coefficient is more significant, the predictions of the model are quantitatively reliable, whereas the agreement with the experimental data becomes weaker when the dependence of D on x becomes milder. Once the agreement for the diffusion coefficient has been assessed we focus on the validity of the Fick-Jacobs approximation for what concerns the MFPT between the origin and an arbitrary position x: Interestingly, the predictions of Eq. (27) match well with the data reported in Ref. [45] in the case of more corrugated channels whereas it is just qualitative for more shallow shapes. Interestingly, this trend is similar to that of the model (based on finite-element numerical simulations) reported in Ref. [45]. Next, we characterize the dependence of the ratio between the maximum and the minimum local diffusion coefficients on the length of the rod. As expected, Fig. 3 shows that in the asymptotic limits, L/h 0 ≪ 1 and L/h 0 ≫ 1, the rods either reduce to a spherical particle of radius l = L or to a very long stiff filament that, due to the confinement is almost always aligned with the channel axis. In both cases the effective diffusion coefficient reduces to D ≃ D 0 . In contrast, for L ≃ h 0 , Fig. 3 shows that increasing the corrugation of the channel increases the range of values of the rod length for which the local diffusion coefficient is sensitive to the local channel section. Clearly, if the aspect ratio L/l approaches unity the diffusion coefficient becomes homogeneous along the channel.

Free energy barrier
Having an explicit formula for the free energy, it allows us to discuss the dependence of the equilibrium (f = 0) free energy barrier defined as the difference between the free energy at the bottleneck and the one at the channel's widest section: For comparison we recall that the free energy difference for an ideal gas depends solely on the geometry of the channel: Our model shows that the dependence of the free energy barrier of rods strongly depends on their length. Indeed, while for shorter rods, L/h 0 0.5, the free energy barrier of the rod increases upon increasing channel corrugation (i.e. increasing ∆A gas ) for longer rods, L/h 0 1, the free energy barrier of the rod decreases upon increasing ∆A gas , as shown in Fig. 4. Additionally, for every geometry of the channel the maximum departure of ∆A from ∆A gas is attained when the length of the rod is comparable to the channel average section, L ≃ h 0 . Interestingly, for rods much longer than the channel section we have that ∆A → 2∆A gas . In particular, very long rods L ≫ h 0 imply that z ≃ 0, i.e. that the center of mass is confined close to the channel axis. Hence, this asymptotic behavior can be understood by expanding Eq. (24) about z = 0. Indeed, for L/h 0 ≫ 1 the integrand of Eq. (24) can be approximated by h(x)/L and hence we obtain that A(x) ∝ 2 ln [h(x)] which eventually leads to ∆A → 2∆A gas .

Transport
Next, we analyze the transport of rods under the action of a constant force. In order to simplify the analysis we assume that the fluid rods are suspended and keep at rest, i.e. there is no advection of rods due to fluid motion. In such a regime, Fig. 5 shows that the dependence of the dimensionless permeability µ on the length of the major axis of the rod, L, is quite complex. Indeed, for flat channels (β∆A gas = 0) the transport is determined solely by the effective diffusion coefficient. In particular, for L ≫ h 0 the rod is almost aligned with the axis of the channel and the effective diffusion coefficient approaches the one of the minor axis (see Fig.3) and, according to Eq.(22) µ → 1 for L ≫ h 0 . In contrast, for β∆A gas = 0 and L ≪ h 0 the rod can freely rotate and it experiences a reduction in the effective diffusion coefficient (see Fig.3). For β∆A gas = 0, the channel is not flat and the dependence of the dimensionless permeability µ on L becomes more involved. Indeed, Fig. 5 shows that for intermediate values of β∆A gas , µ displays a non-monotonous dependence on L/h 0 . This behavior is similar to that observed for polymers confined between corrugated plates (see Ref. [40]). Finally, for larger values of β∆A gas , µ monotonically decreases upon increasing L/h 0 .

B. Charged rods
Next, we analyze the case of charged rods. For simplicity we assume that the charge, q, is localized at the center of the rod and that the channel walls are characterized by a constant charge density σ. In addition, we assume that a dilute monovalent electrolyte is suspended in the fluid phase so that the system is electrically neutral. In this regard, the rod brings an extra charge, whose magnitude is assumed to be much smaller than that of the charge of the double layer such that it does not affect significantly the local electrostatic field 1 . Accordingly, the (unperturbed) electrostatic potential inside the channel reads with φ 0 = σ εκ . The local free energy becomes 1 Such an approximation becomes exact in the limit of vanishing charge of the rod. where the diffusion coefficient reads The dependence of the channel permeability on the electrostatic interaction of the rods with the walls is shown in Fig. 6a. In particular, Fig. 6a shows that the dependence of µ on the interaction potential of the rods with the walls is non-monotonous. For a given charge density on the walls, when the magnitude of the charge is very large (for both signs of the charge) the permeability becomes vanishingly small. Indeed, for such cases the effective free energy barrier to be overcome is very large and hence µ becomes very small. In contrast, upon reducing the magnitude of the charge, µ increases and it reaches its maximum for weakly negatively charged rods, i.e. when rods are attracted by the channel walls. As shown in Fig. 6a the optimal value of the charge depends on the channel corrugation and it decreases upon increasing the channel corrugation. This non-monotonous dependence, and in particular the location of the maximum, can be understood by looking at the equilibrium free energy difference β∆A eq . Indeed, Fig. 6b shows that the maxima observed in Fig. 6a occur for values of q for which β∆A eq = 0 (see Fig. 6b) i.e., when the equilibrium overall free energy barrier drops and the rods have only to overcome smaller local free energy barriers. Fig. 6a shows that the geometry of the channel can be used to tune the sorting of rods depending on their charge. In order to exploit the geometry of the channel and to separate charged rods depending on their length, it is mandatory to find the set of parameters for which the sensitivity of the dimensionless channel permeability on the channel geometry is maximized. As discussed, our model highlights a direct correlation between the equilibrium free energy barrier ∆A eq and the dimensionless channel permeability µ. This implies that the region of the parameter space in which µ is likely to be sensitive to the geometry, ∆A gas , can be guessed by inspecting the equilibrium free energy difference ∆A eq . Fig. 7 shows that for kh 0 1, the isolines β∆A eq = 0 are parallel to the ordinate axis for L h 0 . Hence, all rods with length L h 0 are expected to have relatively similar values of µ. Accordingly, for this set of parameters only rods with length L h 0 can be separated. In contrast, for kh 0 ≃ 0.5 and for β∆A gas ≃ 0.75 the contour line β∆A eq = 0 bends. This can be the signal of a strong dependence of the velocity on the rod size. Fig. 8a shows that for larger values of β∆A gas of and βqφ 0 −0.1, µ is sensitive to L/h 0 only for L h 0 whereas µ shows a plateau for L h 0 , in agreement with our argument based on β∆A eq . In contrast, for smaller values of β∆A gas , Fig. 8b shows that for βqφ 0 −0.1, µ is more sensitive to L for the full range of values explored in Fig. 8b, again in agreement with our argument. Finally we have studied the dependence of the flux of charged rods as a function of the Debye length. Interestingly, Fig. 9 shows that the flux of negatively charged rods (hence attracted by the channel walls) has a nonmonotonous dependence on the dimensionless inverse Debye length κh 0 and it displays a maximum around κh 0 1 for all the geometries of the channel that we have explored. In contrast, the net flow of positively charged rods increases monotonously with the inverse Debye length κh 0 .

C. Conclusions
We have studied the dynamics of charged rods embedded in varying-section channels. Under the assumption of slowly varying channel sections, ∂ x h(x) ≪ 1 we have extended the Fick-Jacobs approximation to the case of charged rods. Our approximation allows us to derive an expression for the local diffusion coefficient, D(x) (see Eq. (17)). We have tested our prediction for the case of neutrally charged rods against experimental and numerical results (see Ref. [45]). Interestingly, our predictions (see Fig. 2) match well with the experimental/numerical data for both the local diffusion coefficient as well as for the Mean First Passage Time. We remark that, while in experiments hydrodynamic coupling between the rods and the channel walls is naturally accounted for, it is not so in our model. Hence, a first result of our analysis is that, in the regime under study, hydrodynamic coupling plays a minor role in the diffusion of rods within corrugated channels. In order to grasp the relevance of the rod geometry on the free energy barrier, ∆A, we have compared it to that of point particles, ∆A gas . Interestingly, we found that the enhancement in the free energy difference, ∆A, is at most twofold (see Fig. 4c). Next, we have characterized the channel permeability, µ, to neutral rods as a function of their length and found that for intermediate values of β∆A eq , µ displays a non-monotonous dependence (see Fig. 5) on the major axis of the rod, L, similarly to what has been observed for confined polymers [40].
Finally, we have characterized the channel permeability to charged rods. In particular, we have focused on the dependence of µ on the rod size and charge. Interestingly, we found that µ displays a non-monotonous dependence on the rod charge and that the value of the charge that maximizes µ is the one for which the overall equilibrium free energy barrier, ∆A eq vanishes. Hence, our model allows us to rationalize the non-monotonous dependence of µ with respect to the charge of the rod. We have exploited this feature to identify those set of parameters for which µ is likely to be more sensitive to the length, L, of the major axis of the rods. Interestingly, this approach has revealed to be reliable in identifying those regimes for which the sensitivity of µ on L is maximized. All in all our results show that µ can be tuned by properly combining the channel geometry and the length and charge of the rod. Moreover, our model allows to rationalize the appearance of non-monotonous dependence of µ on both, the charge of the rod and its length. Our results can be useful for the design of novel micro-and nano-fluidic devices aiming at sorting stiff filaments.