# Geometry controlled dispersion in periodic corrugated channels

Matthieu Mangeat, Thomas Guérin, David S. Dean

► To cite this version:

Matthieu Mangeat, Thomas Guérin, David S. Dean. Geometry controlled dispersion in periodic corrugated channels. EPL - Europhysics Letters, 2017, 118 (4), pp.40004. 10.1209/0295-5075/118/40004. hal-01578776

HAL Id: hal-01578776

https://hal.science/hal-01578776

Submitted on 29 Aug 2017

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
# Geometry controlled dispersion in periodic corrugated channels

M. MANGEAT, T. GUÉRIN and D. S. DEAN

*Laboratoire Ondes et Matière d'Aquitaine (LOMA), CNRS, UMR 5798/Université de Bordeaux F-33400 Talence, France*

PACS 05.40.-a – Fluctuation phenomena, random processes, noise, and Brownian motion

PACS 05.60.Cd – Classical transport

PACS 66.10.cg – Mass diffusion, including self-diffusion, mutual diffusion, tracer diffusion, etc.

**Abstract** – The effective diffusivity $D_e$ of tracer particles diffusing in periodically corrugated axisymmetric two- and three-dimensional channels is studied. The majority of the previous studies of this class of problems are based on perturbative analyses about narrow channels, where the problem can be reduced to an effectively one-dimensional one. Here we show how to analyze this class of problems using a much more general approach which even includes the limit of infinitely wide channels. Using the narrow- and wide-channel asymptotics, we provide a Padé approximant scheme that is able to describe the dispersion properties of a wide class of channels. Furthermore, we systematically identify all the exact asymptotic scaling regimes of $D_e$ and the accompanying physical mechanisms that control dispersion, clarifying the distinction between smooth channels and compartmentalized ones, and identifying the regimes in which $D_e$ can be linked to first passage problems.

**Introduction.** – How fast does a cloud of tracer particles, moving stochastically in a complex heterogeneous medium, disperse? This question naturally appears in a wide range of contexts, including mixing [1–3], sorting [4], contaminant spreading [5] or chemical reactions kinetics [6]. The characterization of dispersion properties, which result from a non-trivial interplay between the geometry of the heterogeneous medium and the transport by forces and/or flows, is an active field of research [1–4,7–11]. At large length and time scales, dispersion is usually characterized by an effective diffusion tensor whose components can be considerably different from typical microscopic diffusivities [12]; canonical examples for increased and decreased diffusivities are given by, respectively, motion in shear hydrodynamic flows (called Taylor dispersion [13]) and in periodic [14] and random [15] potentials.

Here, we consider diffusion of non-interacting particles in channels of non-uniform cross-section, a paradigm for diffusion in confined environments [16,17], arising in contexts as varied as biological cells [18,19], zeolites, porous media, ion channels and microfluidic devices. It is well known that, in the absence of hydrodynamic flow, the effective diffusivity of particles in channels is lower than the microscopic diffusivity. Qualitatively, this can be understood by considering the entropy $S(z)$, which measures

the number of available lateral configurations at fixed longitudinal position $z$: the narrow regions have a reduced entropy and act as entropic barriers, while the wide regions can be viewed as entropic traps, leading to a motion slower than in a uniform channel.

The first quantitative results on diffusion in channels are attributed to Jacobs [20] who derived the first form of the so-called Fick-Jacobs (FJ) approximation. This standard approach, and its various extensions [21–31], are based on a dimensional reduction, and approximate the dynamics of the tracer longitudinal position $z(t)$ by a diffusive dynamics in an entropic potential $\phi(z) \equiv -TS(z)$, possibly with a position-dependent diffusion coefficient $D(z)$. Once the dimensional reduction is carried out, the effective diffusion constant can be computed using exact one-dimensional results [32–34]. Such FJ-like approaches rely however on the assumption that the equilibration dynamics of the lateral position is fast compared to longitudinal motion, which unavoidably leads to a limited range of validity. It has been recognized that the case of abrupt changes of channel radius requires an improvement the one-dimensional description at the cost of employing more sophisticated methods [26,35,36]. A different picture, in principle valid for channels constituted of pores separated by narrow necks, relies on the assumption that the motion is controlled by the first passage events of tracer particles
<figure><img src="image_2.png" /><figcaption>Fig. 1: Schematic of a two-dimensional channel of local width $2R(z)$, or slice of an axisymmetric channel in three dimensions of local radius $R(z)$. The complete channel is formed by the periodic repetition of this motif.</figcaption></figure>

between pores. Calculations of effective diffusivities that rely on first passage time (FPT) arguments have been so far restricted to particular simplified geometries, such as sinusoidal channels [37], septate channels (made of perfectly cylindrical connected cavities) [38,39], or channels formed by overlapping circles [40] or spheres [41]. In general however, the regimes of validity of FJ-like approximations and FPT-approaches are different, and it is therefore difficult to describe the transition between these regimes (except for the sinusoidal channel [37]).

In the present paper, we revisit theoretically the problem of dispersion in two- and three-dimensional axisymmetric channels of arbitrary shape. Our approach uses an exact formula of the effective diffusivity, expressed in terms of an auxiliary function that satisfies a set of partial differential equation at the scale of a single period, which we analyze using singular perturbation analysis and conformal mapping techniques. We systematically identify all the (exact) asymptotic scaling regimes of diffusivity and the accompanying physical mechanisms that control dispersion. In many cases, especially the case of highly corrugated channels, the dispersion coefficient is found to depend on only a few quantities related to the channel geometry rather than on the full details of its shape. We show how the identification of regimes far outside the validity of the one-dimensional effective description can lead to an accurate description of the effective diffusivity for a wide range of parameters, via a Padé approximant. We identify the regimes in which $D_e$ is linked to FPT problems. We also show that, depending on the behavior of the radius near the neck, we can classify channels into smoothly and highly corrugated ones, for which the effective diffusivity displays qualitatively different behaviors.

**Channel geometry and general equations for the effective diffusivity.** – We consider here the problem of the diffusion of an overdamped particle, of microscopic diffusivity $D_0$, in a two- or three-dimensional axisymmetric channel (fig. 1), assumed to be periodic with period $L$. We denote by $z$ the (longitudinal) position in the direction parallel to the channel axis, and we assume that the channel radius $R(z)$ is parametrized as

$$R(z) = a + Hg(z/L), \quad (1)$$

where $a$ is the minimal channel radius, $H$ is the amplitude of variation of the channel radius, and $g$ is a dimensionless periodic function of period 1 which describes the geometrical shape of the channel boundaries, chosen to have a maximal value equal to 1 and a minimal value 0. We define the dimensionless parameters, which we will see determine the various modes of dispersion,

$$\xi \equiv H/a, \quad \varepsilon \equiv a/L. \quad (2)$$

Channels of uniform width thus correspond to $\xi = 0$, while $\xi$ is large for highly corrugated ones. The limit of weakly varying channels thus correspond to $\varepsilon \to 0$ (at fixed $\xi$). Finally we denote by $\Omega$ the unit periodic cell, and we call $V$ its volume.

We aim to characterize the long time effective diffusion coefficient of tracer particles $D_e \equiv \lim_{t \to \infty} [z(t) - z(0)]^2/(2t)$, where the overbar denotes ensemble average. The starting point of our analysis is the following exact expression:

$$D_e = D_0 \left( 1 + \frac{(d-1) \langle f_S R^{d-2} \partial_z R \rangle}{\langle R^{d-1} \rangle} \right), \quad (3)$$

where the notation $\langle w \rangle = \int_0^L dz w(z)/L$ is used for the uniform average over one period for any function $w$, $d$ is the spatial dimension ($d = 2$ or 3), and $D_e$ is expressed in terms of an auxiliary function $f_S(z) \equiv f(r = R(z), z)$, where $f(r, z)$ satisfies

$$\partial_z^2 f + r^{2-d} \partial_r [r^{d-2} \partial_r f] = 0, \quad (4)$$

$$[(\partial_z R)\partial_z f - \partial_r f]_{r=R(z)} = \partial_z R, \quad (5)$$

$$f(r, z + L) = f(r, z); \quad \partial_r f|_{r=0} = 0, \quad (6)$$

where $r$ is the distance to the central axis. These equations (3)–(6) are a particular case of the general description of dispersion in arbitrary periodic systems introduced in refs. [11,42], they are also compatible with the equations of the macrotransport theory of Brenner and Edwards [12]. They express the macroscopic diffusion coefficient $D_e$ as a function of the microscopic structure of the channel, at the scale of one single period. Such a system of partial differential equations can be readily integrated numerically by using standard finite element solvers, leading to the curves presented in figs. 2 and 3 for various channels. $D_e$ is represented as a function of $\varepsilon = a/L$ for different values of the corrugation parameter $\xi$. These curves clearly display two plateaus separated by an intermediate regime; we will now study these asymptotic regimes analytically.

**Slowly varying channels ($\varepsilon \to 0$).** – The first limiting case to consider is that of a slowly varying channel, which here corresponds to the limit $\varepsilon \to 0$, a limit in which the FJ approximation applies since equilibration in the perpendicular direction is much faster than in the longitudinal direction. At leading order, a tracer particle exhibits the effectively one-dimensional dynamics of a Brownian particle $z(t)$ with diffusion coefficient
<figure><img src="image_3.png" /><figcaption>Fig. 2: (Color online) Effective diffusivity for the bidimensional channel of radius $R(z) = a\{3/2 + 0.266[\cos(2\pi z/L) + \sin(6\pi z/L)]\}$ for small (a) and finite (b) values of $\varepsilon$. The channel shape is represented in inset. On both plots, disks represent the numerical values of $D_e/D_0$ obtained by solving eqs. (3)–(6), and continuous lines correspond to the Padé approximant (14). In (a), the first orders of the expansion of $D_e$ in powers of $\varepsilon$, obtained from refs. [23,43], are represented. In (b), we also represent the results obtained by using one-dimensional re-summed formulas for the local diffusivity $D(z)$ proposed by Zwanzig (Zw) [21], Reguera and Rubi (RR) [22] and Kalinay and Percus (KP) [23].</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Fig. 3: (Color online) Effective diffusivity $D_e$ for channels of sinusoidal shape $g(u) = [1 + \cos(2\pi u)]/2$ in two dimensions (a) and three dimensions (b), and ellipsoidal shape $g(u) = \sqrt{1 - 4u^2}$ in two dimensions (c) and three dimensions (d). Disks represent the numerical solution of eqs. (3)–(6), continuous lines correspond to the Padé approximant (14). Dashed lines represent the various asymptotic regimes: FJ expression (7) for $\varepsilon \to 0$, wide channel limit (8) for $\varepsilon \to \infty$ and narrow-escape regime (21) for intermediate $\varepsilon$ (the value $\kappa = 2/\pi$, valid for $H \gg L$, was used).</figcaption></figure>

$D_0$ advected by the potential $\phi(z) = -k_B T \ln(R^{d-1}(z))$. Here, the Lifson-Jackson formula [32] provides an estimate for the effective diffusion coefficient $D_e$:

$$D_e = \frac{D_0}{\langle (1 + \xi g)^{d-1} \rangle \langle (1 + \xi g)^{1-d} \rangle} \equiv D_{\text{FJ}}. \quad (7)$$

This well known expression clearly shows (from Jensen's inequality) that the effective diffusivity $D_e$ is reduced compared to the microscopic diffusion coefficient $D_0$, it is

furthermore independent (at leading order) of the channel period $L$. This estimate can be recovered from the equation for $f$ by a standard perturbation theory in $\varepsilon$ for $d = 2$ [43], and we show in the Supplementary Material Supplementarymaterial.pdf (SM) of this paper how to generalize to $d = 3$.

Several works have attempted to improve this estimate, using various approaches. The most obvious one consists of calculating more terms in the expansion in $\varepsilon$: this has been done by assuming that the dynamics for $z(t)$ can be described by a Markovian one, with a position-dependent local diffusivity $D(z)$. Perturbation expansions for $D(z)$ have been proposed<fnref n="1" /> which have been found to be consistent with the expansion of the macrotransport theory performed up to order $\varepsilon^4$ [43]. However, such series in powers of $\varepsilon$ fail to describe the numerical curve as soon as $\varepsilon$ is not small (see fig. 2), for the obvious reason that at large $\varepsilon$ the curve should reach a plateau instead of being polynomial. The use of Padé approximants is a standard way to enforce a series expansion to have a constant limit at large $\varepsilon$, while retaining precision for small $\varepsilon$: it consists of writing $D_e = \sum_{n=1}^q a_n \varepsilon^n / \sum_{n=1}^p b_n \varepsilon^n$, with $p = q$ in order to ensure a finite limit for large $\varepsilon$, while the coefficients $a_n, b_n$ are chosen to be consistent with the small $\varepsilon$ expansion. We have tried this procedure, but we concluded that it does not lead to accurate results, as the plateau at large $\varepsilon$ is not predicted correctly.

Other approaches [21–23] have considered different choices of $D(z)$, obtained by partial re-summation techniques, and leading to alternative estimates of $D_e$. However it is seen on fig. 2 that none of these re-summations correctly estimate $D_e$ for finite values of $\varepsilon$, and it is also known that they are not consistent with exact small $\varepsilon$ expansion [43]. Therefore, FJ-like approaches are, by construction, not likely to be able to estimate $D_e$ for finite values of $\varepsilon$, which is why we focus on the opposite limit, $\varepsilon \to \infty$ of fast varying channels.

**The limit of wide channels ($\varepsilon \to \infty$).** — In the limit of wide channels, where $a, H \gg L$, the diffusivity at leading order can be deduced as follows. At the time scale $\tau \sim L^2/D_0$, particles at $r < a$ can be considered to diffuse freely in the longitudinal direction, while particles at $r > a$ can be considered as immobile. We can thus estimate the mean square displacement during a time $t$ to be $\overline{z^2(t)} = 2D_0T_c(t)$, where $T_c(t)$ is the time spent in the region $r < a$ up to time $t$. Ergodicity implies that $T_c(t)/t$ is also the ratio of the volume of the region $r < a$ to the total volume of the periodic cell, which leads to

$$D_e = D_0 \frac{a^{d-1}}{\langle (a + Hg)^{d-1} \rangle}. \quad (8)$$

<footnote n="1">¹In ref. [23] an anisotropy of the microscopic diffusion tensor is considered, and the small parameter of the perturbation expansion is the ratio $D_\parallel/D_\perp$; expansions in powers of $\varepsilon$ or of this small parameter are equivalent.</footnote>
This expression is the same as that found in comb-like geometries [44,45] or tubes with dead-end regions [46] in simplified geometries. However, this argument should hold only for infinitely thin dead-end regions, and it does not take into account tracer particles that cross the hypersurface $r = a$, and corresponding corrections to the effective diffusivity are not easy to estimate. In what follows we carry out a quantitative analysis of the exact equations (4) in the large $\epsilon$ limit.

In order to construct the auxiliary function $f$ in the limit of wide channels $\epsilon \to \infty$, it is convenient to use rescaled variables, $\tilde{z} = z/L$ and $\tilde{r} = r/a$, in which case the variation range of variables $\tilde{z}, \tilde{r}$ is independent of $\epsilon$. At leading order in $\epsilon$, the resulting equation (4) for $f(\tilde{r}, \tilde{z})$ becomes $\partial_{\tilde{z}}^2 f = 0$, which, using the boundary conditions (5) and (6) leads to solutions of the form

$$f(\tilde{r}, \tilde{z}) = \tilde{z} \theta(\tilde{r} - 1)L + b(\tilde{r}), \quad (9)$$

where $\theta$ is the Heaviside function and $b(\tilde{r})$ is an undertermined function of $\tilde{r}$. The above solution for $f$ is not satisfactory because it is not continuous at $\tilde{r} = 1$. This is the signal of the presence of a boundary-layer near $r = a$. The size of the boundary layer for the lateral variable $r$ is found by inspection to be $L$: it corresponds to the region in which tracer particles can cross the line $r = a$ in the time $L^2/D_0$ needed by tracer particles to reach a neighboring pore. It is now useful to write $r = a + \eta L$, in which case the equation (4) becomes at leading order in $\epsilon$

$$\partial_{\eta}^{2} f + \partial_{\tilde{z}}^{2} f = 0, \quad (10)$$

$$\begin{cases} \partial_{\tilde{z}} f|_{\tilde{z}=\pm 1/2} = L & (\eta > 0), \\ f(\eta, \tilde{z} = 1/2) = f(\eta, \tilde{z} = -1/2) & (\eta < 0), \end{cases} \quad (11)$$

(note that this equation holds in dimensions 2 and 3). Furthermore, to match with the outer solution (9), $f$ must behave as $f \simeq b(1) + \tilde{z}L$ for $\eta \to \infty$, and $f$ must be constant for $\eta \to -\infty$. This problem can now be handled by the use of complex analysis: we look for a solution $f = \text{Re}(w(Z))$, where $w$ is an analytic function of the complex variable $Z = \tilde{z} + i\eta$. If we make the transformation $Z_1 = ie^{i\pi Z}$, the problem becomes equivalent to the two-dimensional electrostatic problem consisting of finding the potential generated by two perfectly conducting neighboring horizontal plates, being located between $(\pm 1, 0)$ and $(\pm \infty, 0)$, on which opposite values of the potential is imposed. The solution of this problem can be constructed using a Schwarz-Christoffel transform (see SM), and we find

$$f = \text{Re} \left[ \frac{iL}{\pi} \ln \left( 1 + \sqrt{1 + e^{-2\pi i \tilde{z} + 2\pi \eta}} \right) \right] + b(1), \quad (12)$$

where $\text{Re}(\dots)$ represents the real part of a complex number. It can be checked that the above formula satisfies the boundary conditions (11) and matches with the outer solution (9) when one takes $\eta \to \pm\infty$. Inserting this formula

into eq. (3) yields

$$D_e = \lim_{\epsilon \to \infty} D_0 \frac{a^{d-1}}{\langle (a + Hg)^{d-1} \rangle} \left( 1 + \frac{(d-1) \ln 2}{\pi \epsilon} \right), \quad (13)$$

where the $\epsilon^{-1}$ correction comes from the contribution of $f$ in the boundary layer. These corrections, which quantify the contribution to dispersion of the particles that can cross the line separating the blocked region from the regions of free longitudinal move, do not depend on the details of the channel geometry: they are characterized by a universal numerical constant equal to $\ln 2/\pi$.

**An approximant including both narrow- and wide-channel limits.** – At this stage, we can construct a Padé type approximant for $D_e$,

$$D_e = D_{FJ} \frac{1 + a_1 \epsilon + a_2 \epsilon^2 + a_3 \epsilon^3}{1 + b_1 \epsilon + b_2 \epsilon^2 + b_3 \epsilon^3}, \quad (14)$$

where the coefficients $a_i, b_i$ are carefully chosen to ensure that the expression for $D_e$ is exact for both the wide-channel limit $\epsilon \to \infty$ (up to order $\epsilon^{-1}$, using eq. (13)) and the slowly varying channel limit $\epsilon \to 0$ (up to order $\epsilon^4$, for which we used expressions in the literature [43]), see SM. This approximant incorporates effects that cannot be captured by FJ-like approaches, and is found to agree with the numerical curve for almost all values of $\epsilon$ (see figs. 2 and 3). We therefore emphasize that the strength of our approach is that it allows an accurate description of $D_e$ which ranges from narrow to wide channels.

**The FJ approximation for highly corrugated channels.** – We now proceed to simplifying the description of the mechanisms controlling dispersion in the limit of large ratio $\xi = H/a$ of maximum width over minimal aperture. Consider first the large $\xi$ limit of the FJ expression. The result of taking $\xi \to \infty$ in eq. (7) depends on the existence of the integral $\int dz/g^{d-1}(z)$, which may be a divergent one (because $g$ vanishes for some value of $z$). We now assume that the behavior of $R$ near the point of minimal aperture (here taken as the origin of longitudinal axis) is characterized by

$$R(z \to 0) \simeq a + \gamma |z|^{\nu}, \quad (15)$$

where $\gamma$ is a quantity that characterizes the local geometry of the narrowest region of the channel. For example, differentiable channel profiles correspond to $\nu = 2$, in which case $\gamma$ is half the minimal curvature at the neck. If the neck is composed of connected conical portions (so $\nu = 1$), arctan($\gamma$) is half the opening angle of these cones. The assumption (15) is equivalent to

$$g(\tilde{z} \to 0) \simeq A |\tilde{z}|^{\nu}; \quad \gamma = AH/L^{\nu}. \quad (16)$$

If we define $\nu_c(d) = 1/(d-1)$, we see that the integral of $1/g^{d-1}$ is infinite when $\nu > \nu_c$. In this case the dominant contribution in the integral $J \equiv \int_0^1 d\tilde{z}/(1+\xi g(\tilde{z}))^{d-1}$
comes from the values of $z$ close from the points of smallest channel width, so that $J \approx \int_{-\infty}^{\infty} d\tilde{z}/(1+\xi A|\tilde{z}|^\nu)^{d-1}$ (where we can replace the integration bounds by $\pm\infty$ without changing the integration result at leading order). Computing this integral leads to

$$ \frac{D_{\text{FJ}}}{D_0} \approx \frac{\nu \sin(\pi/\nu)(A\xi)^{1/\nu}}{2\pi \xi^{d-1}\langle g^{d-1} \rangle} \left(\frac{\nu}{\nu-1}\right)^{d-2}, \quad (17) $$

which can also be written as

$$ D_{\text{FJ}} \approx L^2/(2T), \quad (18) $$

with

$$ T = \frac{V}{D_0 a^{d-1-1/\nu} \gamma^{1/\nu}} \times \frac{\pi}{2\nu \sin(\pi/\nu)} \left( \frac{2(\nu-1)}{\pi\nu} \right)^{d-2}. \quad (19) $$

We can interpret $T$ as the mean first time to reach the middle of one of the narrow regions, while the other is reflecting. The time $T$ does not depend on the precise geometric details of the channel shape: it depends only on the volume $V$ of a single pore, on the minimal channel radius $a$ and on the parameter $\gamma$ which characterizes the geometry of the channel near the neck. In this regime, the stochastic trajectories of the tracer particles can be viewed as a continuous time random walk, where the particles spend in each pore an average time $T/2$ which measures the rate at which the tracer particles can escape the entropic barriers formed by the narrow regions. Equation (17) is known in the case $\nu = d = 2$ [47,48]. The mean escape time $T$ to an opening at the end of a funnel has recently been calculated using conformal mapping techniques [19,49–51] for $\nu = 2$ and $d = 2, 3$ and coincides with the above formula<fnref n="2" />, it is interesting to see that these mean escape times are also accessible via the FJ approximation. The general formula (19) for $T$ for any exponent $\nu$ is new to the best of our knowledge.

An important remark here is that $D_e$ is controlled by the time to cross the neck regions: as a consequence, eq. (17) holds as soon as the FJ approximation is a correct description of the dynamics in the neck only rather than in the whole channel. The relevant longitudinal length scale $l^*$ in the neck is identified from $a \sim \gamma(l^*)^\nu$, so that $l^* \sim (a/\gamma)^{1/\nu}$; the FJ approximation is valid when $l^* \gg a$, a condition which is less constraining (for $\nu > 1$) than the condition $H \ll L$ which would be required for the FJ approximation to hold in the whole channel.

Thus, if $\nu > \nu_c$, the dispersion in the limit of slowly varying channels is controlled by the geometry at the neck. The situation is completely different in the case $\nu < \nu_c$ for which the large $\xi$ limit of $D_{\text{FJ}}$ reads

$$ D_{\text{FJ}} = \frac{D_0}{\langle g^{d-1} \rangle \langle g^{1-d} \rangle}. \quad (20) $$

In this case, the effective diffusivity depends on the channel’s geometrical shape, but not on any of the parameters $a, L, H$. This is a key difference between channels with sharp necks ($\nu < \nu_c$) or smooth necks ($\nu > \nu_c$): dispersion in sharp neck channels is not controlled by the diffusion at the neck only. Interestingly, the case $\nu = 1$ in $d = 3$ dimensions is included in the regime $\nu > \nu_c$ and corresponds to a regime where the dynamics at the neck controls the transitions between pores and thus the dispersion.

**Intermediate regime of dispersion.** – We finally study the regime that is intermediate between the limits of small and large $\epsilon$. It is seen on fig. 3 that this intermediate regime tends to increase with increasing $\xi$, and also tends to deviate from the predictions of our Padé approximant. This suggests the presence of a different mechanism that controls dispersion. We treated this case by performing a singular perturbation analysis of eq. (4)–(6) in the limit of small pore opening by following closely the approach of refs. [52,53] (see SM for details). We obtain

$$ D_e \approx \frac{L^2 D_0}{V} \times \begin{cases} \frac{2a}{\pi} & (d=3), \\ \frac{\ln(2L\kappa/a)}{2\ln(2L\kappa/a)} & (d=2), \end{cases} \quad (21) $$

where $\kappa$ is a constant that depends on the ratio $H/L$ and on the shape of the boundary; more precisely $\ln \kappa = [R(\mathbf{r}_0, \mathbf{r}_0) + R(\mathbf{r}_1, \mathbf{r}_1) - 2G(\mathbf{r}_0, \mathbf{r}_1)]\pi/2$, where $G$ is the pseudo-Green’s function of the domain (without opening), $R$ is the non-diverging part of this Green’s function and $\mathbf{r}_0, \mathbf{r}_1$ are the positions of the openings. In the limit $H \gg L$, $\kappa$ reaches a constant value deduced from the Green’s function in an infinite strip [54], $\kappa = 2/\pi$. The above formula reveals that in this intermediate regime one can again interpret the stochastic trajectories as continuous time random walks, with a dispersion coefficient satisfying the relation (18), $D_e = L^2/(2T)$. In 3 dimensions, $T$ is, not surprisingly, the mean escape time through a small opening embedded in a flat plane, which does not depend on the initial position of the walker, due to the non-compact feature of space exploration by a Brownian walker in 3D [6,55]. In 2D the situation is slightly different, because Brownian search for an opening is only marginally compact, and mean escape times depend logarithmically on the initial position [6,55]. Comparing eq. (21) with recent calculations of the mean escape time in 2D domains of arbitrary shape [53] reveals that $D_e = L^2/(2T)$, where $T$ is not the global mean first passage time to a pore, but is instead the time to reach a pore, starting from the opposing opening (considered as reflecting). The above formula has been identified for particular geometries such as septate channels in 3D [38] and for channels made of overlapping spheres [41], it has already proposed for the corresponding cases in 2D [39,40] but at leading order only.

We obtained the formula (21) rigorously from (3) in the limit of small pore opening in the case $\nu \le 1$, but one can see from fig. 3 that it is actually gives a good description of $D_e$ in the intermediate regime for large $\xi$ for arbitrary

<footnote n="2">²Note that, for $d = 3$ and $\nu = 2$, eq. (17) is half the result given in refs. [19,50]. It is mentioned in ref. [51] that a correction factor of one half should be added, but misprints in the definition of $R$ and $a$ render difficult the comparison with (17).</footnote>
<figure><img src="image_5.png" /></figure>

Fig. 4: Diagram $(\epsilon, \xi)$ summarizing the asymptotic estimates of $D_e$ and their validity regimes for $\nu > \nu_c$ (a) and $\nu \le \nu_c$ (b), with $\nu_c = 1/(d-1)$. In the Fick-Jacobs (FJ), Narrow-Escape (NE) and Wide-Channel (WC) domains, the equations (7), (21) and (8) are respectively correct. In the limit $\xi \to 0$, the effective diffusivity goes to $D_0$ for both FJ and WC limits. When $\nu > \nu_c$, the FJ expression for large $\xi$ corresponds to a regime where dispersion is controlled by events of narrow escape through a funnel (NEF).

geometries for any channel shape, be it smooth or not. It is therefore not limited to compartmentalized channels. This can be understood by noting that the large $\xi$ limit implies that the boundaries become more and more perpendicular to the channel axis near the channel necks, and one therefore recovers the conditions of the narrow-escape problem at a domain boundary.

We end our study by drawing in fig. 4 qualitative diagrams where the asymptotic expressions for $D_e$ are summarized, together with their validity domains. Each regime corresponds to a different physical mechanism that controls the behavior of the stochastic trajectories and thus dispersion. We stress that our approach, based on the exact expression (3) for $D_e$, enables to obtain all the asymptotic regimes.

**Conclusion.** – Let us now summarize our findings. Here we have revisited the problem of computing the effective diffusivity of tracer particles in corrugated axisymmetric two- and three-dimensional channels. We have classified the channels into two categories: smooth channels, characterized by an exponent $\nu > 1/(d-1)$, for which the FJ dispersion becomes controlled by the crossing of a funnel at the necks, which we computed for any $\nu$, and non-smooth channels, with $\nu < 1/(d-1)$, for which the effective diffusivity in the FJ regime becomes independent of the parameters $H, L, a$ in the strong corrugation limit. We also identified two supplementary regimes, common to all channel geometries: a comb-like regime for wide channels, where we quantified the influence on dispersion of the probability of crossing the frontier between the slow and fast regions, and an intermediate regime controlled by the standard narrow-escape problem. We have also proposed a Padé type approximant for $D_e$, which accurately describes the effective diffusivity for a wide class of parameters between the limits of narrow and wide channels. This study thus provides a refined understanding of how dispersion properties are controlled by the geometry of the channel.

## REFERENCES

[1] LE BORGNE T., DENTZ M. and VILLERMAUX E., *Phys. Rev. Lett.*, **110** (2013) 204501.

[2] DENTZ M., LE BORGNE T., ENGLERT A. and BIJELJIC B., *J. Contam. Hydrol.*, **120** (2011) 1.

[3] BARROS F. P., DENTZ M., KOCH J. and NOWAK W., *Geophys. Res. Lett.*, **39** (2012).

[4] BERNATE J. A. and DRAZER G., *Phys. Rev. Lett.*, **108** (2012) 214501.

[5] BRUSSEAU M. L., *Rev. Geophys.*, **32** (1994) 285.

[6] CONDAMIN S., BÉNICHOU O., TEJEDOR V., VOITURIEZ R. and KLAFTER J., *Nature*, **450** (2007) 77.

[7] LEITMANN S. and FRANOSCH T., *Phys. Rev. Lett.*, **118** (2017) 018001.

[8] AMINIAN M., BERNARDI F., CAMASSA R., HARRIS D. M. and MCLAUGHLIN R. M., *Science* (2016) 0532.

[9] HAYNES P. and VANNESTE J., *J. Fluid Mech.*, **745** (2014) 321.

[10] TZELLA A. and VANNESTE J., *Phys. Rev. Lett.*, **117** (2016) 114501.

[11] GUÉRIN T. and DEAN D. S., *Phys. Rev. Lett.*, **115** (2015) 020601.

[12] BRENNER H., *Macrotransport Processes* (Elsevier) 2013.

[13] TAYLOR G., *Proc. R. Soc. Lond. A*, **219** (1953) 186.

[14] VAN KAMPEN N., *Stochastic Processes in Physics and Chemistry*, third edition (North-Holland, Amsterdam) 2007.

[15] DEAN D. S., DRUMMOND I. and HORGAN R., *J. Stat. Mech.: Theory Exp.* (2007) P07013.

[16] BURADA P. S., HÄNGGI P., MARCHESON F., SCHMID G. and TALKNER P., *ChemPhysChem*, **10** (2009) 45.

[17] MALGARETTI P., PAGONABARRAGA I. and RUBI M., *Front. Phys.*, **1** (2013) 21.

[18] BRESSLOFF P. C. and NEWBY J. M., *Rev. Mod. Phys.*, **85** (2013) 135.

[19] HOLCMAN D. and SCHUSS Z., *Rep. Prog. Phys.*, **76** (2013) 074601.

[20] JACOBS M., *Diffusion Processes* (Springer, New-York) 1967.

[21] ZWANZIG R., *J Phys. Chem.*, **96** (1992) 3926.

[22] REGUERA D. and RUBI J., *Phys. Rev. E*, **64** (2001) 061106.

[23] KALINAY P. and PERCUS J., *Phys. Rev. E*, **74** (2006) 041203.

[24] KALINAY P. and PERCUS J., *Phys. Rev. E*, **72** (2005) 061203.

[25] KALINAY P. and PERCUS J., *J. Chem. Phys.*, **122** (2005) 204701.

[26] KALINAY P. and PERCUS J. K., *Phys. Rev. E*, **82** (2010) 031143.

[27] MARTENS S., SCHMID G., SCHIMANSKY-GEIER L. and HÄNGGI P., *Phys. Rev. E*, **83** (2011) 051135.

[28] BRADLEY R. M., *Phys. Rev. E*, **80** (2009) 061142.

[29] BEREZHKOVSKII A. and SZABO A., *J. Chem. Phys.*, **135** (2011) 074108.

[30] DAGDUG L. and PINEDA I., *J. Chem. Phys.*, **137** (2012) 024107.

[31] VALDES C. V. and GUZMAN R. H., *Phys. Rev. E*, **90** (2014) 052141.

[32] LIFSON S. and JACKSON J. L., *J. Chem. Phys.*, **36** (1962) 2410.
[33] REIMANN P., VAN DEN BROECK C., LINKE H. HÄNGGI P., RUBI J. M. and PÉREZ-MADRID A., *Phys. Rev. Lett.*, **87** (2001) 010602.

[34] REGUERA D., SCHMID G., BURADA P. S., RUBI J. M., REIMANN P. and HÄNGGI P., *Phys. Rev. Lett.*, **96** (2006) 130603.

[35] BEREZHKOVSKII A. M., BARZYKIN A. V. and ZITSERMAN V. Y., *J. Chem. Phys.*, **131** (2009) 224110.

[36] ANTIPOV A. E., BARZYKIN A. V., BEREZHKOVSKII A. M., MAKHNOVSKII Y. A., ZITSERMAN V. Y. and ALDOSHIN S. M., *Phys. Rev. E*, **88** (2013) 054101.

[37] BOSI L., GHOSH P. K. and MARCHESONI F., *J. Chem. Phys.*, **137** (2012) 174110.

[38] BORROMEO M. and MARCHESONI F., *Chem. Phys.*, **375** (2010) 536.

[39] MARCHESONI F., *J. Chem. Phys.*, **132** (2010) 166101.

[40] PINEDA I., VAZQUEZ M.-V., BEREZHKOVSKII A. M. and DAGDUG L., *J. Chem. Phys.*, **135** (2011) 224101.

[41] BEREZHKOVSKII A. M., ZITSERMAN V. Y. and SHVARTS-MAN S. Y., *J. Chem. Phys.*, **118** (2003) 7146.

[42] GUÉRIN T. and DEAN D. S., *Phys. Rev. E*, **92** (2015) 062103.

[43] DORFMAN K. D. and YARIV E., *J. Chem. Phys.*, **141** (2014) 044118.

[44] BEREZHKOVSKII A. M., DAGDUG L. and BEZRUKOV S. M., *J. Chem. Phys.*, **142** (2015) 134101.

[45] DEAN D. S. and JANSONS K. M., *J. Stat. Phys.*, **70** (1993) 1313.

[46] DAGDUG L., BEREZHKOVSKII A. M., MAKHNOVSKII Y. A. and ZITSERMAN V. Y., *J. Chem. Phys.*, **127** (2007) 224712.

[47] DAGDUG L., VAZQUEZ M.-V., BEREZHKOVSKII A. M., ZITSERMAN V. Y. and BEZRUKOV S. M., *J. Chem. Phys.*, **136** (2012) 204106.

[48] CRANK J., *The Mathematics of Diffusion* (Oxford University Press) 1979.

[49] HOLCMAN D., HOZE N. and SCHUSS Z., *Phys. Rev. E*, **84** (2011) 021906; *Phys. Rev. E*, **85** (2012) 039903 (erratum).

[50] HOLCMAN D. and SCHUSS Z., *Multiscale Model. Simul.*, **10** (2012) 1204.

[51] GUERRIER C. and HOLCMAN D., *Eur. Phys. J. ST*, **223** (2014) 3273.

[52] WARD M. J. and KELLER J. B., *SIAM J. Appl. Math.*, **53** (1993) 770.

[53] PILLAY S., WARD M. J., PEIRCE A. and KOLOKOLNIKOV T., *Multiscale Model. Simul.*, **8** (2010) 803.

[54] BARTON G., *Elements of Green's Functions and Propagation* (Clarendon Press, Oxford) 1989.

[55] BÉNICHOU O. and VOITURIEZ R., *Phys. Rev. Lett.*, **100** (2008) 168105.