Biased diffusion in confined media: Test of the Fick-Jacobs approximation and validity criteria

We study biased, diffusive transport of Brownian particles through narrow, spatially periodic structures in which the motion is constrained in lateral directions. The problem is analyzed under the perspective of the Fick-Jacobs equation which accounts for the effect of the lateral confinement by introducing an entropic barrier in a one dimensional diffusion. The validity of this approximation, being based on the assumption of an instantaneous equilibration of the particle distribution in the cross-section of the structure, is analyzed by comparing the different time scales that characterize the problem. A validity criterion is established in terms of the shape of the structure and of the applied force. It is analytically corroborated and verified by numerical simulations that the critical value of the force up to which this description holds true scales as the square of the periodicity of the structure. The criterion can be visualized by means of a diagram representing the regions where the Fick-Jacobs description becomes inaccurate in terms of the scaled force versus the periodicity of the structure.

tion, being based on the assumption of an instantaneous equilibration of the particle distribution in the cross-section of the structure, is analyzed by comparing the different time scales that characterize the problem. A validity criterion is established in terms of the shape of the structure and of the applied force. It is analytically corroborated and verified by numerical simulations that the critical value of the force up to which this description holds true scales as the square of the periodicity of the structure. The criterion can be visualized by means of a diagram representing the regions where the Fick-Jacobs description becomes inaccurate in terms of the scaled force versus the periodicity of the structure.

I. INTRODUCTION
In many transport phenomena such as those taking place in biological cells, ion channels, nano-porous materials and microfluidic devices etched with grooves and chambers, Brownian particles, instead of diffusing freely in the host liquid phase, undergo a constrained motion.
The uneven shape of these structures regulates the transport of particles yielding important effects exhibiting peculiar properties. The results have implications in processes such as catalysis, osmosis and particle separation [1,2,3,4,5,6,7,8,9] and for the noise-induced transport in periodic potential landscapes that lack reflection symmetry (ratchet systems) [10,11,12] or ratchet transport mechanisms which are based on asymmetric geometries, termed "entropic" ratchet devices [10,11,12,13]. For example, it has been found that the separation of DNA fragments in narrow channels [14,15,16] is largely influenced by their shape. The translocation of structured polynucleotides through nanopores also allows one to determine their sequence and structure [17,18,19].
The motion of the particles through these quasi-one-dimensional structures can in principle be analyzed by means of the standard protocol of solving the Smoluchowski equation with the appropriate boundary conditions imposed. Whereas this method has been very successful when the boundaries of the system possess a regular shape, the challenge to solve the boundary value problem in the case of uneven boundaries represents typically a very difficult task. A way to circumvent this difficulty consists in coarsening the description by reducing the dimensionality of the system, keeping only the main direction of transport, but taking into account the irregular nature of these boundaries by means of an entropic potential. The resulting kinetic equation for the probability distribution, the Fick-Jacobs (F-J) equation, is similar in form to the Smoluchowski equation, but now contains an entropic term. The entropic nature of this term leads to a genuine dynamics which is very different from that observed when the potential has an energetic origin [20]. It has been shown that the F-J equation can provide a very accurate description of entropic transport in 2D and 3D channels of varying cross-section [20,21].
However, the derivation of the F-J equation entails a tacit approximation: The particle distribution in the transverse direction is assumed to equilibrate much faster than in the main (unconstrained) direction of transport. This equilibration justifies the coarsening of the description leading in turn to a simplification of the dynamics, but raises the question about its validity when an external force is applied. To establish the validity criterion of a F-J description for such biased diffusion in confined media is, due to the ubiquity of this situation, a subject of primary importance.
Our objective with this work is to investigate in greater detail the F-J approximation for biased diffusion and to set up a corresponding criterion describing its regime of validity. We will analyze the biased movement of Brownian particles in 2D and 3D periodic channels of varying cross section and formulate different criteria for the validity of such a F-J description.
The paper is organized in the following way: In Sec. II, we describe the physical situation and introduce the model defined through the corresponding Langevin and Fokker-Planck equations. In Sec. III, we introduce the F-J approach for the unbiased situation and extend it to the driven case. Sec. IV is devoted to establish a criteria for the validity of the F-J approximation derived by comparing the different characteristic time scales of the process.
In Sec. V, the accuracy of the F-J description is tested against numerical simulations for a 2D periodic channel, and the conditions of validity are summarized in a diagram in terms of the scaled force versus the periodicity of the structure. In Sec. VI, we provide further explanations on when and why the equilibration assumption fails and the F-J approach leads to useable results. Finally, in Sec. VII we present our main conclusions.

II. DIFFUSION IN CONFINED SYSTEMS
In typical transport processes through pores or channels (like the one depicted in Fig.   1), the motion of the suspended particles is induced by application of an external potential V ( r) resulting in a force F . In general, the dynamics of the suspended Brownian particles is governed by Langevin's equation: where r is the two or three dimensional position vector of a particle of mass m, η is its friction coefficient, k B the Boltzmann constant, T the temperature, and a dot over the quantity refers to time derivative. The thermal fluctuations due to the coupling of the particle with the environment are modeled by a zero-mean Gaussian white noise ξ(t), obeying the fluctuationdissipation relation ξ i (t) ξ j (t ′ ) = 2 δ ij δ(t − t ′ ) for i, j = x, y, z. In the over-damped case, The half-width ω is a periodic function of x with periodicity L.
i.e. when m¨ r(t) << η˙ r(t), the inertia term in Eq. (1) can safely neglected and the Langevin equation describing the dynamics of a Brownian particle within the channel reads: In addition to Eq. (2), the full problem is set up by imposing reflecting boundary conditions at the channel walls.
The corresponding Fokker-Plank equation for the time evolution of the probability distribution P ( r, t) takes the form [22,23]: where J( r, t) is the probability current: and denotes the diffusion coefficient of the suspended particles.
Due to the impenetrability of the channel walls, the normal component of the probability current J( r, t) vanishes at the boundaries. If n denotes the vector perpendicular to the channel walls, the reflecting boundary conditions read: In this paper we focus on the case of a symmetric 2D channel where the force is constant and directed along the axis, cf. Fig. 1. The half-width of the 2D channel is given by a periodic function ω(x), i.e. ω(x + L) = ω(x) for all x. In this case, the boundary condition at y = ±ω(x). For an arbitrary form of ω(x), the boundary value problem defined through Eqs. (3), (4) and (7) is very difficult to solve. Despite the inherent complexity of this problem an approximate solution can be found by introducing an effective one-dimensional description where geometric constraints and bottlenecks are considered as entropic barriers [9,20,21,24,25,26].

III. THE FICK-JACOBS APPROXIMATION
In the absence of an external force, i.e. when F = 0, it was shown [24,25] that the dynamics of particles in confined structures (such as that of Fig. 1) can be described by the (the prime refers to the first derivative, i.e. ω ′ (x) = dω(x)/dx), and it has been shown that the introduction of an effective x-dependent diffusion coefficient can considerably improve the accuracy of the kinetic equation, thus extending its validity to more winding structures [21,25,26]. The expression: where α = 1/3, 1/2 for two and three dimensions, respectively, has been shown to accurately account for the curvature effects [21].
In the presence of a constant force F along the direction of the channel the F-J equation can be recast into the expression [20,21] −F x denotes the energy contribution and S = k B ln h(x) the entropy contribution. For a symmetric channel with periodicity L, the free energy assumes the form of a periodic tilted potential.
Note that Eq. (10) may typically describe the time-evolution of a 1D particle distribution within an energy landscape. In the present context, however, due to the reduction of the geometric confinements in 2D or 3D space into one dimension, we end up with an entropic contribution to the free energy. In absence of a fixed bias F = 0, we deal with a pure entropic situation with the free energy reading (10) Recently, we have shown that the dynamics of a confined Brownian particle, in presence of an applied bias, can accurately be described by means of equation (10) [20]. In the presence of very strong applied bias, and for more winding structures, however, the F-J equation becomes inaccurate. In the present work we present further numerical and analytical results and will set up tailored criteria under which the F-J approximation assumes good validity.

A. Particle current
One of the key quantities in transport through quasi-one-dimensional structures is the study of the average particle current ẋ , or equivalently the nonlinear mobility, which is defined as the ratio of the average particle current and the applied force F . For the average particle current we derive an expression which is similar to the Stratonovich formula for the current in titled periodic energy landscapes [27,28,29,30], but with a spatial diffusion coefficient. A detailed derivation of this expression is given in the Appendix A, cf. Eq. (A12).
Hence, we obtain the nonlinear mobility for a 2D channel with a shape defined by ω(x): where we have made use of Eq. (9) and the relations D 0 = k B T /η and β = 1/k B T . Substituting z ′ = x ′ /L and z = x/L, the nonlinear mobility scaled with the friction coefficient η can be expressed, for the case of constant forcing, in terms of a single, dimensionless scaling parameter [20]: Therefore, Eq. (12) for where Eq. (14) can be transformed into Eq. (6) in Ref. [20] by interchanging the order of integration.
The asymptotic values of the nonlinear mobility for the cases f → 0 and f → ∞, can be evaluated to read: and where is the average over a spatial period, given an arbitrary periodic function g(x).
Similarly, the time scale associated to diffusion in the axial direction is The third time scale is defined through the characteristic time associated to the drift (ballistic motion) over a distance ∆x given by ∆x where we have used the scaling factor f = F L/k B T . Rearranging the previous expression, we obtain In order to have a good equilibration in the transverse direction, the characteristic time scale associated to diffusion in this direction has to be much smaller than the other two time scales. Therefore, in the absence of an external force, equilibration in the transverse direction occurs if τ y /τ x ≪ 1. This results in the condition: which constitutes the validity criterion of the F-J approach such as put forward by Zwanzig [25].
In presence of a force along the axis, equilibration in the transverse direction demands that the condition τy τ drif t ≪ 1 also holds. Consequently, where in the second step we have replaced the characteristic distances ∆y by the width 2ω(x), and ∆x by L.
A general estimate of the criteria that has to be satisfied is that max( τy τx , τy τ drift ) ≪ 1. An even stronger criteria in order for the F-J description to hold in presence of a constant force can be put forward by considering the sum of the two ratios, cf. Eqs. (22) and (23), i.e.: Eq. (24) provides a quite stringent criteria that indicates when the F-J description of a system is expected to be valid. Note also that this is a local criterion, i.e. for a given channel, there would be regions (associated to drastic changes in the shape of the channel, i.e ω ′ (x) 2 ≫ 1 ) where equilibration in the transverse direction is not feasible, whereas in others is fulfilled. It is then more convenient to work with a global criterion of validity rather than with a local one. One way of getting that global condition is by averaging the local criterion over the period L of the channel, yielding one of our main results: In order to get an explicit estimate of the dependence of the maximum force value on the periodicity L of the channel, we define a critical force value f c , for which the inequality (25) becomes an equality, i.e. ω ′ (x) 2 + 2fc L 2 ω(x) 2 = 1. Then the critical force value reads: which indicates that the critical force scales asymptotically as L 2 , if we fix the overall shape of the channel and change only its periodicity.
Eq. (26) provides an estimate of the minimum forcing beyond which the F-J description is expected to fail in providing an accurate description of the dynamics of system. The quantitative value of the critical force will obviously depend on the level of accuracy sought.
What is really important is how it depends (or scales) with the relevant parameters of the problem.

V. NUMERICAL SIMULATIONS FOR A 2D CHANNEL
In order to check the consistency of the criteria proposed and the validity of the F-J description in the presence of a force, we will compare the analytical result for the scaled nonlinear mobility obtained in Eq. (14) with the simulation results of the overdamped Langevin dynamics in Eq. (2) for the 2D periodic channel sketched in Fig. 1. We remark here that the extension of our scheme to 3D with a rotational symmetry along the transport axis is possible as well. This will consume more computation time, but the overall findings remain qualitatively robust. This feature has been verified with a few numerical tests.
The shape of the channel is described by In order to test the accuracy of the F-J description, we have evaluated the behavior of the nonlinear mobility as a function of the scaled force f , for different values of the periodicity L.
In Fig. 2 For our example of a 2D channel whose shape is defined by Eq. (27), the validity criterion given by Eq. (26) simplifies to thus predicting that the critical value of the force scales as L 2 .
This prediction has been verified by the simulations. Fig. 3 shows the value of the critical force for a tolerance of 1% as a function of the periodicity L. The critical value of the force depends quadratically on the periodicity L 2 , as predicted.
In Fig. 4 we illustrate, for the considered two-dimensional channel, the regions where the F-J approximation is accurate when compared with the simulations. We depict the dependence on the periodicity L of the maximal critical scaled force, obtained by comparing numerical results with the analytic solution for the particle current, for different required relative errors. This diagram shows the regions in parameter space in terms of f and L for which an accurate solution is provided. Thus, it is possible to provide an accurate result by using the analytic solution over a wide range of the scaled parameter and the periodicity. the scaled force f . As the force increases, we can clearly see how the particles are not homogeneously distributed in y-direction, evidencing the failure of the equilibration assumption.
This effect is specially dramatic in Fig. 5(lower right panel), where the force is so strong that the particles do not fully explore the available space in y-direction.
A more detailed analysis could be provided by checking the normalized steady-state probability distribution in the transverse direction at a given x-positions, i.e.
P st x (y) := In Fig. 6 only feel the presence of the boundaries when they are close to the bottlenecks. Hence, in the limit of very large force values, the influence of the entropic barriers practically disappears.
In this limit, the correction in the diffusion coefficient leading to a spatial dependency, i.e.
In the case of a periodic diffusion coefficient, D(x + L) = D(x), and a tilted periodic free energy, A(x + L) = A(x) − F L, it is convenient to define the reduced probability density and the corresponding current, J(x, t) = n J(nL + x, t) , n ∈ Z .
By definition these functions are periodic with periodicity L,P (x + L, t) =P (x, t) and J(x + L, t) =Ĵ(x, t). MoreoverP (x, t) andĴ(x, t) enter a continuity equation (A2) and P (x, t) is normalized on any interval (x, x + L), provided that P (x, t) is normalized, e.g. . (A10) The general relation between the stationary probability current and the steady state par- which implies ẋ =Ĵ L. Then, the particle current reads: Remarkably, this expression for the particle current is equivalent to the expression obtained via the mean-first-passage-time approach presented in Ref. [20].