Hostname: page-component-797576ffbb-xg4rj Total loading time: 0 Render date: 2023-12-05T08:39:58.285Z Has data issue: false Feature Flags: { "corePageComponentGetUserInfoFromSharedSession": true, "coreDisableEcommerce": false, "useRatesEcommerce": true } hasContentIssue false

Interactions enhance dispersion in fluctuating channels via emergent flows

Published online by Cambridge University Press:  27 September 2023

Y. Wang
Courant Institute of Mathematical Sciences, New York University, New York, NY, USA
D.S. Dean
University of Bordeaux, CNRS, LOMA, UMR 5798, F-33400 Talence, France Team MONC, INRIA Bordeaux Sud Ouest, CNRS UMR 5251, Bordeaux INP, University of Bordeaux, F-33400 Talence, France
S. Marbach*
Courant Institute of Mathematical Sciences, New York University, New York, NY, USA CNRS, Sorbonne Université, Physicochimie des Electrolytes et Nanosystèmes Interfaciaux, F-75005 Paris, France
R. Zakine*
Courant Institute of Mathematical Sciences, New York University, New York, NY, USA Chair of Econophysics and Complex Systems, École Polytechnique, 91128 Palaiseau Cedex, France LadHyX UMR CNRS 7646, École Polytechnique, 91128 Palaiseau Cedex, France
Email addresses for correspondence: [email protected], [email protected]
Email addresses for correspondence: [email protected], [email protected]


Understanding particle motion in narrow channels can guide progress in numerous applications, from filtration to vascular transport. Thermal or active fluctuations of fluid-filled channel walls can slow down or increase the dispersion of tracer particles via entropic trapping in the wall bulges or hydrodynamic flows induced by wall fluctuations, respectively. Previous studies concentrated primarily on the case of a single Brownian tracer. Here, we address what happens when there is a large ensemble of interacting Brownian tracers – a common situation in applications. Introducing repulsive interactions between tracer particles, while ignoring the presence of a background fluid, leads to an effective flow field. This flow field enhances tracer dispersion, a phenomenon reminiscent of that seen for single tracers in incompressible background fluid. We characterise the dispersion by the long-time diffusion coefficient of tracers numerically and analytically with a mean-field density functional analysis. We find a surprising effect where an increased particle density enhances the diffusion coefficient, challenging the notion that crowding effects tend to reduce diffusion. Here, inter-particle interactions push particles closer to the fluctuating channel walls. Interactions between the fluctuating wall and the now-nearby particles then drive particle mixing. Our mechanism is sufficiently general that we expect it to apply to various systems. In addition, our perturbation theory quantifies dispersion in generic advection–diffusion systems.

JFM Papers
Creative Commons
Creative Common License - CCCreative Common License - BYCreative Common License - SA
This is an Open Access article, distributed under the terms of the Creative Commons Attribution-ShareAlike licence (, which permits re-use, distribution, and reproduction in any medium, provided the same Creative Commons licence is used to distribute the re-used or adapted article and the original article is properly cited.
© The Author(s), 2023. Published by Cambridge University Press.

1. Introduction

Fluid transport in confined channels, and generally in porous structures, is relevant for a broad range of biological and industrial applications, from nutrient transport in vascular networks and microorganisms in soils (Bhattacharjee & Datta Reference Bénichou, Illien, Oshanin, Sarracino and Voituriez2019; Tomkins, Hughes & Morris Reference Tomkins, Hughes and Morris2021), to improved filtration of liquids or gases, including modern desalination techniques and oil recovery (Werber, Osuji & Elimelech Reference Werber, Osuji and Elimelech2016; Marbach & Bocquet Reference Marbach and Bocquet2019). However, narrow channels present a number of challenging features to model, e.g. the predominance of surface effects, the importance of spatiotemporal fluctuations, as well as specific electrostatic response (Kavokine, Netz & Bocquet Reference Kavokine, Netz and Bocquet2021). Significant progress in the last decade has improved our understanding of transport features in such porous environments, and we briefly review below the effects of (i) spatial and (ii) temporal fluctuations of the confining environment, as well as (iii) crowding effects due to interactions.

First, purely spatial corrugations of the confining channel reduce the long-time diffusion coefficient of isolated particles. In fact, random crossings of channel constrictions are rare events that impede overall transport: the constrictions form effective entropic barriers (Zwanzig Reference Zwanzig1992). This effect is well captured by the so-called Fick–Jacobs formalism (Jacobs Reference Jacobs1935; Reguera & Rubí Reference Reguera and Rubí2001; Kalinay & Percus Reference Kalinay and Percus2006; Rubi Reference Rubi2019), which reduces the problem to an effectively one-dimensional one. The Fick–Jacobs formalism has been extended recently to arbitrary channel geometries (Dagdug, García-Chung & Chacón-Acosta Reference Dagdug, García-Chung and Chacón-Acosta2016; Chávez, Chacón-Acosta & Dagdug Reference Chávez, Chacón-Acosta and Dagdug2018). Recent work suggests that the approach is also adapted to study fluid flow through biological membranes (Arango-Restrepo & Rubi Reference Arango-Restrepo and Rubi2020). Such entropic contributions induce significant corrections to transport in microfluidic (Yang et al. Reference Yang, Liu, Li, Marchesoni, Hänggi and Zhang2017) or biological (Rubí et al. Reference Rubí, Lervik, Bedeaux and Kjelstrup2017) channels.

Second, and of importance in several applications, confining channels are often not static but fluctuate in time, due to either thermal agitation or an external forcing. Molecular dynamics simulations found, early on, enhanced diffusion of gas in microporous materials if the thermal vibrations of the material are accounted for (Leroy, Rousseau & Fuchs Reference Leroy, Rousseau and Fuchs2004; Haldoupis et al. Reference Haldoupis, Watanabe, Nair and Sholl2012), and more recently, enhanced water diffusion in solid state pores via phonon-fluid coupling (Ma et al. Reference Ma, Grey, Shen, Urbakh, Wu, Liu, Liu and Zheng2015, Reference Ma, Tocci, Michaelides and Aeppli2016; Cao, Wang & Ma Reference Cao, Wang and Ma2019; Noh & Aluru Reference Noh and Aluru2021). Recent theoretical work suggests that enhanced diffusion is universally due to longitudinally fluctuating fluid flows, driven by fluctuations of the channel walls, that convect the tracer particles (Marbach, Dean & Bocquet Reference Marbach, Dean and Bocquet2018). The mechanism is thus reminiscent of Taylor–Aris dispersion, where the long-time diffusion coefficient is enhanced by a cross-sectionally inhomogeneous fluid flow profile (in the initially studied case of Poiseuille flow) (Taylor Reference Taylor1953; Aris Reference Aris1956; Hydon & Pedley Reference Hydon and Pedley1993). When the characteristic time scale of wall fluctuations is smaller (larger) than the time scale to diffuse across typical constrictions, diffusion is enhanced (decreased), a criterion that was verified numerically (Chakrabarti & Saintillan Reference Chakrabarti and Saintillan2020) and in experiments in fluctuating porous matrices (Sarfati, Calderon & Schwartz Reference Sarfati, Calderon and Schwartz2021). Such enhancement of dispersion properties is relevant in a number of biological contexts, such as in blood vessels (Masri, Puelz & Riviere Reference Masri, Puelz and Riviere2021), slime mould vasculature (Marbach & Alim Reference Marbach and Alim2019), the gut (Codutti, Cremer & Alim Reference Codutti, Cremer and Alim2022), and near molecular motors (Evans, Krause & Feringa Reference Evans, Krause and Feringa2021), and are of general relevance in plants (Tomkins et al. Reference Tomkins, Hughes and Morris2021).

Third, in addition to the fluctuating confinement, tracer particles are often not isolated and interact with other particles or molecules that are also diffusing in the medium. The current paradigm is that such crowded environments, in an open domain, tend to slow down diffusion at equilibrium (Lekkerkerker & Dhont Reference Lekkerkerker and Dhont1984; Lowen & Szamel Reference Lowen and Szamel1993; Dean & Lefèvre Reference Dean and Lefèvre2004). However, subtle effects may arise if the crowded medium is driven out of equilibrium (Zaccone et al. Reference Zaccone, Dorsaz, Piazza, De Michele, Morbidelli and Foffi2011; Lappala, Zaccone & Terentjev Reference Lappala, Zaccone and Terentjev2013). For example, the diffusion coefficient of a tracer driven by an external force may be enhanced at low density and high forcing (Bénichou et al. Reference Bénichou, Illien, Oshanin and Voituriez2013, Reference Bertini, De Sole, Gabrielli, Jona-Lasinio and Landim2018; Démery, Bénichou & Jacquin Reference Démery, Bénichou and Jacquin2014; Illien et al. Reference Illien, Bénichou, Oshanin, Sarracino and Voituriez2018). Typically, a trade-off is observed between, on the one hand, the increased diffusion induced by faster exploration of space thanks to the driving force, and on the other hand, the decreased diffusion due to spatial constrictions induced by the confining media (Illien et al. Reference Illien, Bénichou, Oshanin, Sarracino and Voituriez2018). Such effects may be exploited for active microrheology within spatially corrugated channels (Puertas, Malgaretti & Pagonabarraga Reference Puertas, Malgaretti and Pagonabarraga2018; Malgaretti, Puertas & Pagonabarraga Reference Malgaretti, Puertas and Pagonabarraga2022).

Given its obvious importance, especially in biological systems, the coupling between the effects of fluctuating channels and inter-particle interactions has received surprisingly little attention. Since temporal channel fluctuations increase transport coefficients, and since inter-particle interactions, or crowding effects, generally decrease diffusion, it is natural to ask if one can predict which effect dominates. Even though recent numerical work showed that diffusion enhancement could be obtained with increased particle density in microporous matrices (Obliger et al. Reference Obliger, Valdenaire, Ulm, Pellenq and Leyssale2019), sometimes even exhibiting a maximum (Pireddu et al. Reference Pireddu, Pazzona, Demontis and Załuska-Kotur2019) at a certain density, it was also noticed that the effect can depend strongly on the precise details of the system understudy (Obliger et al. Reference Obliger, Bousige, Coasne and Leyssale2023). These results thus call for a general theoretical investigation.

In this paper, we investigate the motion of diffusing particles with repulsive interactions in a confined and fluctuating channel (see figure 1a), which is essentially a spatially periodic profile moving with constant velocity $v_{wall}$. We perform numerical simulations to quantify the effective long-time self-diffusion coefficient $D_{eff}$ and the effective long-time drift $V_{eff}$ of particles. We explain their behaviour for a broad range of interaction strengths between particles and fluctuating channel speeds $v_{wall} = \omega _0/k_0$, where $k_0$ and $\omega _0$ are the wavenumber and frequency of the wall fluctuations. Using perturbation theory, we obtain simple analytical predictions for ${V_{eff}}$ and $D_{eff}$ that are in excellent agreement with simulations.

Figure 1. Set-up to study the transport of tracers in an interacting system with fluctuating boundaries. (a) Tracer particles (yellow) perform a random walk within this wiggling environment, here represented by a top, moving wall (red arrows). (b) We consider pairwise, soft, repulsive interactions between particles, and we vary the interaction strength $\alpha$ and the particle number density $\rho _0$ to investigate more or less crowded environments.

The paper is organised as follows. To start, in § 2, we consider a simple ideal gas in the fluctuating channel. We find that the long-time diffusion coefficient $D_{eff}$ exhibits a maximum with respect to the wall phase velocity $v_{wall}$, corresponding to maximal enhancement of diffusion. This result should be contrasted with that of the diffusion of tracer particles in incompressible fluids that exhibits a monotonic increase with the wall phase velocity $v_{wall}$ (Marbach et al. Reference Marbach, Dean and Bocquet2018). The non-monotonic behaviour seen in the case studied here originates from the interplay between diffusion and advection due to the wall that increases long-time diffusion only when the diffusive time scale and the advection time scale are comparable. In § 3, we then consider soft-core interactions between particles (see figure 1b). We find that $D_{eff}$ can be enhanced further with increasing repulsive interactions. Here, repulsive interactions play a role in generating a more uniform distribution of particles in the channel, even in the vicinity of rapidly fluctuating bulges. Eventually, increased wall–particle collisions, caused by now-nearby particles, enhance dispersion. This behaviour is reminiscent of the mechanism of enhanced tracer diffusion in a fluctuating channel filled with an incompressible fluid (Marbach et al. Reference Marbach, Dean and Bocquet2018), and indeed in § 4, we show that remarkably, analytically and numerically, transport coefficients converge to a universal incompressible fluid regime in the limit of strong repulsive interactions. Finally, in § 5 we discuss how these mechanisms and techniques may be used further to investigate more complex situations of transport in fluctuating confined environments.

2. Transport of ideal (isolated) tracers in a fluctuating channel

2.1. Simulation results

We start by exploring the motion of tracers in a fluctuating channel, where the environment is not crowded, i.e. where tracer particles are sufficiently far away from each other that we can consider that they do not interact. The tracer particles perform a random walk with diffusion coefficient $D_0$. We refer to this case as the ‘ideal gas’ regime. Throughout this study, we will assume that the impact of the fluctuating channel boundaries on the velocity field of the supporting fluid is negligible. This is the case of particles embedded in a highly compressible fluid or gas, for which the mean free path is much smaller than any relevant length scale in the system. Examples of such compressible systems in biology include airflow in the pulsating alveoli in the lungs (Sznitman et al. Reference Sznitman, Heimsch, Heimsch, Rusch and Rösgen2007; Dong, Yang & Zhu Reference Dong, Yang and Zhu2022). Alternatively, this situation can also be achieved considering that the boundaries are not implemented by hard walls but rather by potential barriers, such as (fluctuating) electrostatic potentials for charged tracers. The impact of fluctuating boundaries on the fluid's velocity field has already been characterised in the Stokes flow regime (Marbach et al. Reference Marbach, Dean and Bocquet2018).

We simulate the motion of these non-interacting particles in a fluctuating channel, described through a fluctuating wall at height $h(x,t) = H + h_0 \cos ( k_0 x - \omega _0 t )$. The shape of the fluctuating interface is chosen to be sinusoidal, as a generic interface profile can be decomposed in terms of plane waves. The presence of the walls is incorporated via a soft boundary potential acting on the particles. We track the motion of particles over long times, and evaluate their effective long-time self-diffusion coefficient $D_{eff}$ and mean drift $V_{eff}$ along the main channel axis $x$ (see Appendix A for details).

Previous work on incompressible fluids has established that a relevant parameter to analyse the system is given by the Péclet number characterising the fluctuations,

(2.1)\begin{equation} Pe = \frac{\omega_0}{D_0 k_0^2} = \frac{\tau_{diff}}{\tau_{wall}}, \end{equation}

which compares the time scale to diffuse across the length of a channel corrugation $\tau _{diff} = 1/D_0 k_0^2$ to the period of the channel fluctuations $\tau _{wall} = 1/\omega _0$ (Marbach et al. Reference Marbach, Dean and Bocquet2018). In simulations, we therefore fix the typical channel corrugation length $L = 2 {\rm \pi}/ k_0$ and vary the wall phase velocity $v_{wall} = \omega _0/k_0$. All parameters are expressed in terms of a time unit $\tau _0$ and a length unit $\ell _0$, which are arbitrary. We present the results (yellow dots) in figures 2(a,b) for $D_{eff}$ and ${V_{eff}}$ with increasing $Pe$.

Figure 2. Transport properties of non-interacting tracers in a periodic driven channel; comparison between a compressible (ideal gas) and an incompressible supporting fluid. (a) Long-time longitudinal diffusion coefficient and (b) drift of particles, with $Pe = \omega _0/D_0 k_0^2 =v_{wall}L/(2{\rm \pi} D_0)$. In (a) and (b), error bars corresponding to one standard over 10 independent runs are smaller than point size, except for (b) small $Pe$, since $V_{eff} \simeq v_{wall} \sim 0.01 \ell _0/\tau _0$ is small compared to the noise level. Ideal gas theory corresponds to (2.18) and (2.19), and incompressible fluid to (2.24) and (2.25). Stationary density profiles of particles for $v_{wall}=0.1 \ell _0/\tau _0$ (or $Pe \simeq 3$) for (c) the ideal gas and (d) an incompressible fluid, for $v_{wall}=0.5 \ell _0/\tau _0$ ($Pe \simeq 16$). Blue arrows show the velocity field in the channel, represented in the referential of the moving wall, with length proportional to magnitude (arbitrary scale). The colour scale for the density profiles is shared in (c) and (d), and yellow (purple) regions indicate regions of high (low) density. Other numerical parameters are $L=200\ell _0$, $H=12 \ell _0$, $h_0=3\ell _0$ and $D_0=1 \ell _0^2/\tau _0$.

Interestingly, the long-time diffusion coefficient $D_{eff}$ exhibits a non-monotonic dependence on the Péclet number $Pe$ (yellow dots in figure 2a), which we interpret in the following way. At low Péclet numbers $Pe \ll 1$, particles move much faster than the wall, $\tau _{wall} \gg \tau _{diff}$, hence the wall appears frozen. Particles therefore spend time exploring the bulges before escaping the constrictions, and their diffusion coefficient is consequently decreased, $D_{eff} \leq D_0$. This effect is well captured by the Fick–Jacobs approximation, where transport is reduced due to constrictions acting as effective entropy barriers (Jacobs Reference Jacobs1935; Zwanzig Reference Zwanzig1992; Reguera & Rubí Reference Reguera and Rubí2001; Kalinay & Percus Reference Kalinay and Percus2006; Burada et al. Reference Burada, Schmid, Reguera, Rubi and Hänggi2007; Mangeat, Guérin & Dean Reference Mangeat, Guérin and Dean2017; Marbach et al. Reference Marbach, Dean and Bocquet2018; Rubi Reference Rubi2019). At intermediate Péclet numbers $Pe \gtrsim 1$, $\tau _{wall} \simeq \tau _{diff}$ and the moving boundaries enhance particle motion: particles collide with the moving corrugated walls, which increases their diffusion coefficient overall, $D_{eff} \geq D_0$. At high Péclet numbers $Pe \gtrsim 10$, $\tau _{wall} \ll \tau _{diff}$ and the wall moves so fast that particles no longer have time to enter the bulges; they therefore behave as though they were not seeing any corrugation, the effective diffusion coefficient is unchanged and equals the bare diffusion, $D_{eff} = D_0$. In other words, from the point of view of the particles, the channel is effectively flat with height given by its minimum height.

To further understand the behaviour at high Péclet numbers, we plot the average density distribution $\rho (x,y)$ of particles within the channel (the average being over noise), in the reference frame where the channel profile is stationary, for a large value of $v_{wall}$ (figure 2c). We observe that particles accumulate at the constriction. Indeed, in this reference frame, the channel constriction acts as a bottleneck for transport. This can be seen as an inverse Bernoulli effect: in the channel's frame of reference, the particle flux $\rho (x) h(x) v_{wall}$ is necessarily constant, where $\rho (x)$ is the cross-sectional averaged density of particles. Hence the density (and not the velocity, as in the Bernoulli effect) is largest where the channel is constricted, $\rho (x) \propto 1/h(x)$. Eventually, since the particles explore only low vertical coordinates $y \lesssim H - h_0$, they no longer collide with the moving wall, hence their diffusion coefficient is unchanged.

The effective particle drift ${V_{eff}}$ decreases monotonically with the Péclet number (yellow dots in figure 2(b), noting that the vertical axis is the normalised drift relative to the wall phase velocity, ${V_{eff}} /v_{wall}$). At small Péclet numbers, the density of particles is uniform in the channel, therefore the particles within the bulge are carried along in the same direction as the wall phase velocity. At higher Péclet numbers, the fraction of particles within the bulge decreases, as seen in figure 2(c). Since all particles accumulate within the bottleneck, they are no longer pushed by the moving wall, hence ${V_{eff}} /v_{wall} \rightarrow 0$.

2.2. Analytic theory to account for transport in fluctuating channels

To account for this broad behaviour, we build a general analytic theory that reproduces these effects. With the goal of making a pedagogical introduction to our perturbation method, which closely follows that which was described only briefly in Marbach et al. (Reference Marbach, Dean and Bocquet2018), we devote this subsection to a detailed explanation.

2.2.1. Constitutive equations

Brownian tracer particles evolve in a fluctuating environment, confined in the $y$ direction between $y = 0$ and $y = h(x,t)$, and infinitely long in the $x$ direction. Compared to the simulation, we here make no assumption on the form of $h(x,t)$, which may describe thermal motion of the wall (determined e.g. through a Hamiltonian characterising the flexibility of the interface) or active motion (driven by an external source of energy, as in our simulations where fluctuations are imposed). The probability density function $\rho _{tr}(x,y,t)$ to find a tracer particle at position $(x,y)$ at time $t$ obeys the Fokker–Planck equation

(2.2)\begin{equation} \frac{\partial{\rho_{tr}(x,y,t)}}{\partial t}+ \boldsymbol{ \nabla} \boldsymbol{\cdot} \boldsymbol{j} (x,y,t) = 0, \end{equation}

with $\boldsymbol {j}(x,y,t) \equiv - D_0\,\boldsymbol { \nabla } \rho _{tr}(x,y,t) + \boldsymbol {u} (x,y,t)\,\rho _{tr}(x,y,t)$, and where we have assumed that the tracer's diffusion coefficient $D_0$ is uniform in space and that the tracer is advected by the field $\boldsymbol {u} (x,y,t)$, which can have both potential and non-potential components. In the case of ideal non-interacting particles, the underlying supporting fluid does not move due to the moving interface, and there are no interactions, $\boldsymbol {u} (x,y,t) \equiv 0$. For now, we keep $\boldsymbol {u}(x,y,t)$ arbitrary as this will be useful in the incompressible and interacting regimes.

2.2.2. Boundary conditions

We impose no flux boundary conditions at both surfaces. This means that the projection of the current on the direction normal to the lower surface boundary is zero, and that similarly, the projection of the current on the direction normal to the upper boundary in the frame moving with speed $\boldsymbol {U}_\mathrm {wall}$, where the surface is static, is zero. The boundary conditions thus read

(2.3)\begin{gather} j_y(x,0,t) = 0, \end{gather}
(2.4)\begin{gather}\boldsymbol{j} (x,h(x,t),t) \boldsymbol{\cdot} \boldsymbol{n} =- \rho (x,h(x,t),t)\,\boldsymbol{U}_{wall} \boldsymbol{\cdot} \boldsymbol{n}. \end{gather}

Here, $\boldsymbol {n}$ is the outward normal to the interface, and $\boldsymbol {U}_{wall} = (0, \partial _t h)$ is the velocity of the wall, considering that the wall atoms move vertically about their average position, similarly to phonon modes that propagate on material structures (Chaikin, Lubensky & Witten Reference Chaikin, Lubensky and Witten1995) or peristaltic motion of vasculature (Alim et al. Reference Alim, Amselem, Peaudecerf, Brenner and Pringle2013; Marbach & Alim Reference Marbach and Alim2019). We derive these boundary conditions from first principles in § 1 of the supplementary material available at

2.2.3. Simplified longitudinal equation for the probability distribution

To make progress on the long-time behaviour of (2.2), we place ourselves within the lubrication approximation. This means that we consider the typical corrugation length $L$ to be much bigger than the average channel height, $\langle h(x,t) \rangle = H$, $H \ll L$, itself much bigger than the amplitude of the channel fluctuations, i.e. $\sqrt {\langle (h - H)^2 \rangle }\propto h_0 \ll H$. Within the lubrication approximation, the outward normal to the interface in (2.4) simplifies to $\boldsymbol {n} \simeq (\partial _x h,1)$. In this framework, the particle density relaxes much faster on the vertical direction $y$ than on the longitudinal direction $x$. At low Péclet number, this should notably yield a vertically uniform particle density. Figure 2(c) shows that the density profile is indeed independent of $y$ in our simulations. Also, in Appendix B, we show that even though the lubrication approximation should fail for e.g. larger channel heights $H$, the results derived below remain surprisingly robust.

We therefore look for an approximate evolution equation on the probability distribution integrated vertically, or marginal distribution, $p_{tr}(x,t) = \int _0^{h(x,t)} \rho _{tr}(x,y,t) \, \text {d} y$. Taking the time derivative of $p_{tr}(x,t)$ and using (2.2) yields

(2.5) \begin{align} \frac{\partial{p_{tr}(x,t)}}{\partial t} &= \frac{\partial{h(x,t)}}{\partial t}\,\rho_{tr}(x,h(x,t),t)- \int_0^{h(x,t)} \left( \frac{\partial{j_y(x,y,t)}}{\partial y} + \frac{\partial{j_x(x,y,t)}}{\partial x} \right) \text{d} y \nonumber\\ &= \frac{\partial{h(x,t)}}{\partial t}\,\rho_{tr}(x,h(x,t),t) - j_y(x,h(x,t),t) + j_y(x,0,t) \nonumber\\ &\quad + \frac{\partial{h(x,t)}}{\partial x}\,j_x(x,h(x,t),t) - \frac{\partial{}}{\partial x} \left(\int_0^{h(x,t)} j_x(x,y,t) \,\text{d} y\right), \end{align}

where we used the simple and useful relation

(2.6)\begin{equation} \frac{\partial{}}{\partial x} \left( \int_0^{h(x,t)} f(x,y,t) \,\text{d} y\right) = \frac{\partial{h(x,t)}}{\partial x}\,f(x,h(x,t)) + \int_0^{h(x,t)} \frac{\partial{f(x,y,t)}}{\partial x} \,\text{d} y \end{equation}

for any function $f$. Equation (2.5) can be simplified remarkably by using the no-flux boundary conditions (2.4) and (2.3), which lead to all the surface terms vanishing, and we find

(2.7)\begin{equation} \frac{\partial{p_{tr}(x,t)}}{\partial t} =- \frac{\partial{}}{\partial x} \left(\int_0^{h(x,t)} j_x(x,y,t) \,\text{d} y\right). \end{equation}

We look for a closed equation on $p_{tr}(x,t)$. Writing

(2.8)\begin{equation} j_x(x,y,t) =- D_0\,\frac{\partial{\rho_{tr}(x,y,t)}}{\partial x} + u_x\, \rho_{tr}(x,y,t) \end{equation}

explicitly, and using (2.6) again, we find

(2.9)\begin{align} \frac{\partial{p_{tr}(x,t)}}{\partial t}& =- \frac{\partial{}}{\partial x} \left(- D_0\,\frac{\partial p_{tr}(x,t)}{\partial x} + D_0\,\frac{\partial{h(x,t)}}{\partial x}\, \rho_{tr}(x,h(x,t),t)\right.\nonumber\\ &\quad \left.{}+ \int_0^{h(x,t)}u_x\,\rho_{tr}(x,y,t) \,\text{d} y \right), \end{align}

where the last terms still do not depend explicitly on the marginal $p_{tr}$. We will therefore make the common first-order approximation that since the probability distribution profile is nearly uniform vertically, we may assume $\rho _{tr}(x,y,t) \sim p_{tr}(x,t)/h(x,t)$. We then obtain the following Fokker–Planck equation on the marginal distribution, at lowest order in $\rho _{tr}(x,y,t) - p_{tr}(x,t)/h(x,t)$:

(2.10)\begin{equation} \frac{\partial{p_{tr}}}{\partial t} =- \frac{\partial{}}{\partial x} \left(-D_0\,\frac{\partial{p_{tr}}}{\partial x} + D_0 p_{tr}\,\frac{\partial{\ln h}}{\partial x} + \overline{u_x} p_{tr} \right), \end{equation}

where $\overline {u_x} = ({1}/{h}) \int _0^{h(x,t)} u_x \,\text {d} y$ is the vertically averaged longitudinal drift. This final step can be made more rigorous using a centre manifold expansion; see Mercer & Roberts (Reference Mercer and Roberts1990, Reference Mercer and Roberts1994) and Marbach & Alim (Reference Marbach and Alim2019). Equation (2.10) can alternatively be obtained by starting with the ansatz $\rho _{tr}(x,y,t) \sim p_{tr}(x,t)/h(x,t)$ at the level of (2.2), and making use of the boundary conditions. However, the derivation that we present here has the advantage of being approximation-free until (2.9).

Equation (2.10) clearly simplifies the initial problem, and it is sufficient to study the long-time behaviour of $p_{tr}$ to obtain the long-time diffusion coefficient $D_{eff}$ and drift ${V_{eff}}$. When $\overline {u_x} \equiv 0$, (2.10) corresponds exactly with the Fick–Jacobs equation (Jacobs Reference Jacobs1935; Zwanzig Reference Zwanzig1992; Reguera & Rubí Reference Reguera and Rubí2001), which describes the motion of (non-interacting) particles in spatially varying but time independent channels $y \equiv h(x)$. Interestingly, the Fick–Jacobs equation is therefore also valid for fluctuating channels in time, $h(x,t)$, regardless of the functional shape of $h(x,t)$, the only assumption being the lubrication approximation. We note that (2.10) is also consistent with the case of a background incompressible fluid, when $\overline {u_x}$ is the cross-sectionally averaged fluid velocity (Eq. (1) of Marbach et al. Reference Marbach, Dean and Bocquet2018).

2.2.4. Perturbation theory to obtain the long-time transport behaviour

Our goal is now to obtain the long-time transport behaviour of particles from the simplified, longitudinal evolution equation (2.10). We seek an approximate long-time equation for the marginal distribution $p_{eff}(x,t) = p_{tr}(x, t \rightarrow \infty )$ of the form

(2.11)\begin{equation} \frac{\partial{p_{eff}(x,t)}}{\partial t} = D_{eff}\,\frac{\partial^2 p_{eff}(x,t)}{\partial x^2} - V_{eff}\,\frac{\partial{p_{eff}(x,t)}}{\partial x} + \delta(x)\,\delta(t), \end{equation}

such that we can naturally read off the long-time diffusion $D_{eff}$ and drift $V_{eff}$, and where we assumed, without loss of generality, that the particle was initially at the centre of the domain.

To keep the calculations general, we rewrite (2.10) as a Fokker–Planck equation

(2.12)\begin{equation} \frac{\partial{p_{tr}(x,t)}}{\partial t} =- \frac{\partial{}}{\partial x} \left( - D_0\,\frac{\partial{p_{tr}(x,t)}}{\partial x} + v(x,t)\,p_{tr}(x,t) \right) + \delta(x)\,\delta(t), \end{equation}

where $v(x,t)$ is a general advection coefficient (for the ideal gas case, $v = D_0 \partial _x h / h$). Here, we consider that $v(x,t) = O(\varepsilon )$ is a fluctuating perturbation. Its average over realisations of the noise vanishes, $\langle v(x,t)\rangle = 0$, and its fluctuations are described in Fourier space by a spectrum $S(k,\omega )$ as

(2.13)\begin{equation} \langle \tilde{v}(k,\omega)\,\tilde{v}(k',\omega') \rangle = S(k,\omega)\,(2{\rm \pi})^2\,\delta(k + k')\,\delta(\omega + \omega'), \end{equation}

where $k$ and $\omega$ are the wavenumber and frequency, respectively, and the Fourier transform $\tilde {v}(k,\omega )$ of $v(x,t)$ is defined in (C1). Performing a perturbation development to solve (2.12) on the small parameter $\varepsilon$, we find (see Appendix C)

(2.14)\begin{gather} D_{eff} = D_0 \left(1- \varepsilon^2 \int \frac{\text{d} k \,\text{d} \omega}{(2{\rm \pi})^2}\,\frac{k^2( D_0^2 k^4 - 3 \omega^2)}{( D_0^2 k^4+ \omega^2 )^2}\,S(k,\omega)\right), \end{gather}
(2.15)\begin{gather}V_{eff} = \varepsilon^2 \int \frac{\text{d} k \,\text{d} \omega}{(2{\rm \pi}^2)}\,\frac{\omega}{k} \frac{k^2}{D_0^2k^4 + \omega^2}\,S(k,\omega). \end{gather}

Equations (2.14) and (2.15) are two of the main results of this work. They predict the long-time transport properties of particles within a fluctuating channel with arbitrary fluctuating local drift $v(x,t)$. The results can be applied regardless of the nature of the fluctuations, be they thermal or non-equilibrium, and regardless of the strength of the interactions between particles, and between particles and the supporting fluid. In essence, we generalise the results of Marbach et al. (Reference Marbach, Dean and Bocquet2018), which were valid only for an incompressible supporting fluid.

2.3. Applications of the theory to the periodic channel

2.3.1. Ideal gas

To use the above analytic framework to describe our simulations, we now take the case of the periodic travelling wave $h(x,t) = H + h_0 \cos ( 2{\rm \pi} (x - v_{wall} t ) /L ) = H +$ $h_0 \cos ( k_0 x - \omega _0 t)$ on the surface. The Fourier transform of $h - H$ is $\tilde {h}(k,\omega ) = (h_0/2)$ $(\delta (k+k_0)\,\delta (\omega - \omega _0) + \delta (k-k_0)\,\delta (\omega + \omega _0) )$. The local drift $v$ in the ideal gas case is at lowest order in $\varepsilon = h_0/H$,

(2.16)\begin{equation} v(x,t) = D_0\,\frac{1}{H}\,\frac{\partial{h}}{\partial x},\quad \text{and in Fourier space}\quad \tilde{v}(k,\omega) = {\rm i} D_0 k\,\frac{\tilde{h(k,\omega)}}{H}. \end{equation}

This yields a spectrum for the local drift as

(2.17)\begin{equation} S(k,\omega) = \frac{D_0^2 k^2}{(2{\rm \pi})^2}\,\frac{h_0^2}{4 H^2}\,( \delta(k+k_0)\,\delta(\omega - \omega_0) + \delta(k-k_0)\,\delta(\omega + \omega_0) ). \end{equation}

Plugging the expression for the spectrum in (2.14) and (2.15) yields the long-time transport coefficients in the ideal gas case:

(2.18)\begin{gather} \frac{D^{\textit{ideal gas}}_{eff}}{D_0} = 1 + \frac{h_0^2}{2 H^2}\,\frac{3Pe^2-1}{(Pe^2 + 1)^2}, \end{gather}
(2.19)\begin{gather}\frac{V^{\textit{ideal gas}}_{eff}}{v_{wall}} = \frac{h_0^2}{2 H^2}\,\frac{1}{1 + Pe^2}, \end{gather}

where we recall the expression of the Péclet number $Pe = \omega _0/D_0 k_0^2$.

We present plots of (2.18) and (2.19) as lines in figures 2(a,b), and find excellent agreement with simulations. This agreement is robust over a wide range of physical parameters – see figure 9. Analytically, we recover for $Pe \rightarrow 0$ (corresponding to fixed channels, $v_{wall} = 0$) the well-known entropic trapping result where $D_{eff} = D_0 ( 1- h_0^2/2H^2)$ (Zwanzig Reference Zwanzig1992; Jacobs Reference Jacobs1935; Reguera & Rubí Reference Reguera and Rubí2001; Marbach et al. Reference Marbach, Dean and Bocquet2018; Rubi Reference Rubi2019). The analytic computation recovers the non-monotonic behaviour of the diffusion coefficient for intermediate $Pe$, and confirms that $D_{eff} = D_0$ at high $Pe\to \infty$. Finally, the amplitude of the correction for both the diffusion and the drift scales as $h_0^2/H^2$, which can be interpreted naturally and confirms the collision mechanism in bulges; indeed, the strength of the fluctuations scales as $h_0/H$ but only a fraction $h_0/H$ of particles lie within the bulge and take part in the enhancement of the diffusion coefficient or in the mean drift.

2.3.2. Comparison with transport in incompressible fluid

We now relate our results for the ideal gas to transport within incompressible fluids (Marbach et al. Reference Marbach, Dean and Bocquet2018). In that case, when the channel walls fluctuate, because the supporting fluid is incompressible, they induce fluctuations in the fluid's velocity field. Particles are thus advected by fluid flow. The channel walls here are perfectly slipping walls, which corresponds to the smooth boundary conditions considered for the gas particles, which have no specific lateral friction with the walls. The flow field $(u_x,u_y)$ can be derived in the low-Reynolds-number limit, which is the relevant limit to consider since we envision applications in microfluidics to nanofluidics. We find

(2.20)\begin{gather} u_x(x,y,t) = U_0(x,t), \end{gather}
(2.21)\begin{gather}u_y(x,y,t) = \frac{y}{h} \left( \frac{\partial{h(x,t)}}{\partial t} + U_0\,\frac{\partial{h(x,t)}}{\partial x} \right), \end{gather}

where $U_0(x,t)$ is the average fluid flow. We calculate $U_0(x,t)$ assuming peristaltic flow, i.e. that the flow is driven purely by channel fluctuations and that there is no mean pressure-driven flow (Marbach & Alim Reference Marbach and Alim2019; Chakrabarti & Saintillan Reference Chakrabarti and Saintillan2020). The average pressure force on the fluid has to vanish, and we find ${U_0(x,t) \simeq v_{wall} (h_0/H)}$ ${\cos ( 2{\rm \pi} (x - v_{wall} t )/L )}$ at lowest order in $h_0/H$ (see § 2 of the supplementary material). The flow field is presented in figure 2(d) as blue arrows. As the channel moves towards the right-hand side, fluid mass in the right-hand side bulge has to swell out, consistently with the flow lines. Similarly, the bulge on the left-hand side opens up, allowing fluid flow to come in.

The advection term now has two contributions, one coming from the spatial inhomogeneities, and one coming from advection by fluid flow

(2.22)\begin{equation} v(x,t) = D_0\,\frac{1}{H}\,\frac{\partial{h(x,t)}}{\partial x} + U_0(x,t),\end{equation}

or $\tilde {v}(k,\omega ) = ( {\rm i} D_0 k + \omega /k )\,\tilde {h}(k,\omega )/H$ in Fourier space. The spectrum of the fluctuating drift is thus

(2.23)\begin{equation} S(k,\omega) = \frac{ D_0^2 k^2 + \omega^2/k^2 }{(2{\rm \pi})^2}\,\frac{h_0^2}{4 H^2}\,( \delta(k+k_0)\,\delta(\omega - \omega_0) + \delta(k-k_0)\,\delta(\omega + \omega_0) ). \end{equation}

We then obtain

(2.24)\begin{gather} \frac{D_{eff}^{\textit{inc. fluid}}}{D_0} = 1 + \frac{h_0^2}{2H^2}\,\frac{3 Pe^2 -1}{Pe^2 + 1}, \end{gather}
(2.25)\begin{gather}\frac{V_{eff}^{\textit{inc. fluid}}}{v_{wall}} = \frac{h_0^2}{2H^2}. \end{gather}

We perform numerical simulations where particles are advected by the flow field defined by (2.20) and (2.21). We present the numerical results as blue dots and the analytical results as blue lines in figures 2(a,b). We find perfect agreement between simulation and theory, confirming the analytical approach of Marbach et al. (Reference Marbach, Dean and Bocquet2018).

In contrast with the ideal gas, when particles are surrounded by an incompressible fluid, $D_{eff}$ increases monotonically until it reaches a plateau at large $Pe$, and ${V_{eff}}$ stays constant regardless of $Pe$. In fact, at any value of $Pe$, and especially at high $Pe$, the density distribution of particles within the channel is uniform, as can be seen in figure 2(d). Therefore, even at high $Pe$, the population of particles lying within the bulges is pushed by the moving boundaries and increases both $D_{eff}$ and ${V_{eff}}$.

In this case it is natural to ask why the density distribution remains homogeneous along the channel. Looking closely at the velocity field (figure 2d), we see that particles are carried away from the bottleneck by the flow field, and into the bulges. The flow field therefore facilitates recirculation of accumulated particles. This naturally raises the question of how the supporting fluid's compressibility changes the transport properties of particles within fluctuating channels.

3. Interactions increase diffusion and drift in fluctuating channels

To investigate the impact of the compressibility of the supporting fluid, we introduce interactions between the particles (see figure 1b), and tune the interaction strength to vary the effective compressibility of the system.

3.1. Pairwise interactions and compressibility

We simulate the dynamics of interacting particles within a simple periodic fluctuating channel $h(x,t) = H + h_0\cos (k_0 x - \omega _0t)$, as in § 2. We use a pairwise interaction potential between particles,

(3.1)\begin{equation} \mathcal{V}_{int}(r) =\begin{cases} \alpha (r - d_c)^2 & \text{if } r < d_c,\\ 0 & \text{if } r \geq d_c, \end{cases} \end{equation}

where $r$ is the inter-particle distance, $d_c$ is a critical distance characterising the radius of the interaction (see figure 1b), and $\alpha$ is a spring constant that characterises the strength of the interaction. Note that $\alpha >0$ corresponds to repulsive interactions, while $\alpha <0$ corresponds to attractive interactions. We use a soft-core potential (instead of a hard-core potential as in e.g. Bénichou et al. Reference Bénichou, Illien, Oshanin and Voituriez2013; Suárez, Hoyuelos & Mártin Reference Suárez, Hoyuelos and Mártin2015) as it facilitates numerical integration over long time scales, which is necessary to obtain reliable statistics on the diffusion coefficient. Soft-core potentials are used commonly and also simplify analytic treatments (Pàmies, Cacciuto & Frenkel Reference Pàmies, Cacciuto and Frenkel2009; Démery et al. Reference Démery, Bénichou and Jacquin2014; Antonov, Ryabov & Maass Reference Antonov, Ryabov and Maass2021). We will show later that the numerical results are well reproduced by our theory, whose predictions are robust to strong changes of the interaction potential. The mix of interacting particles and the surrounding (compressible) fluid forms a fluid of interacting particles, and we study its properties below.

The choice of interactions is well adapted to tune the compressibility $\chi ^{\star }_T$ of the fluid of interacting particles, and hence probe different compressibility regimes, in between the ideal gas and incompressible fluid cases. Indeed, $\chi ^{\star }_T$ is related to the structure factor $S_f(k)$ of the gas at zero wavelength (Hansen & McDonald Reference Hansen and McDonald2013):

(3.2)\begin{equation} \chi^{{\star}}_{{T}}(\rho_0,\alpha) = \chi_T^{id}\lim_{k \rightarrow 0} S_f(k), \end{equation}

where $\chi _T^{id} = 1/(\rho _0 k_BT)$ is the compressibility of the ideal gas, which is infinite when $\rho _0 \rightarrow 0$. In the so-called and broadly used random phase approximation, one can relate the structure factor $S_f(k)$ to the interaction potential as

(3.3)\begin{equation} S_f(k) = \frac{1}{1 - \rho_0 \left( - \dfrac{1}{k_BT}\,\tilde{\mathcal{V}}_{int}(k)\right)}. \end{equation}

We calculate

(3.4)\begin{equation} \lim_{k\rightarrow 0} \tilde{\mathcal{V}}_{int}(k) = \iint \mathcal{V}_{int}(x,y) \,\text{d} x\, \text{d} y = \frac{\rm \pi}{6}\,\alpha d_c^4. \end{equation}

We therefore obtain the compressibility of the fluid of interacting particles as

(3.5)\begin{equation} \chi^{{\star}}_T(\rho_0,\alpha) =\frac{ \chi_T^{id}}{1 + \dfrac{\rm \pi}{6}\,\dfrac{\rho_0 \alpha d_c^4}{k_B T}} \equiv \frac{ \chi_T^{id}}{1 + \dfrac{E_0}{k_B T}}, \end{equation}

where we introduced $E_0 = ({\rm \pi} /6) \alpha d_c^4\rho _0$, which can be understood as the energy contribution of interactions contained in a typical volume $d_c^d$, with $d$ the dimension of space (here, $d=2$). For $\alpha <0$, one has $E_0<0$, and the fluid of particles ends up with a higher compressibility than the one of the ideal gas. In what follows, we will focus mainly on repulsive interactions ($\alpha >0$), but we keep in mind that our analytic computation does not assume the sign of $\alpha$. For $\alpha >0$, the compressibility $\chi ^{\star }_T(\rho _0,\alpha )$ decreases with the gas density $\rho _0$ and with increasing interaction strength. Therefore, either of the parameters $\rho _0$ and $\alpha$ is a good candidate to probe the intermediate regime between the ideal gas and the incompressible fluid for which $\chi _T^\star \to 0$. In the next subsection, we explore varying values of the density $\rho _0$, and we will find similar results when varying $\alpha$ (reported in Appendix D).

Finally, it is known that the long-time self-diffusion coefficient of interacting particles, in the bulk, i.e. in an open domain, $\overline{D}_0(\rho _0)$, is decreased in a crowded environment in equilibrium, compared to the infinitely dilute case (Démery et al. Reference Démery, Bénichou and Jacquin2014) (although it may be non-monotonic at high densities as the energy landscape is flattened when the local density is large for soft interactions; see figure 7). We therefore first perform simulations within fixed flat walls ($h_0 = 0$, $v_{wall} = 0$), where the system is in equilibrium at temperature $T$. The confinement plays no particular role in that case because our channels are wide enough compared to the typical size of a particle ($d_c \ll H$). We measure the diffusion coefficients $\overline{D}_0(\rho _0)$ for various particle densities in the channel $\rho _0$ and interaction strengths $\alpha$, and our results agree well with existing theories (Démery et al. Reference Démery, Bénichou and Jacquin2014; see Appendix A.4). We can thus now measure changes in the self-diffusion coefficient when there are fluctuations relative to this bulk value $\overline{D}_0(\rho _0)$.

3.2. Increasing interactions enhances diffusion and drift

We perform simulations of interacting particles in fluctuating channels ($h_0 = H/4$, varying $v_{wall} > 0$). We present our numerical results for the long-time self-diffusion coefficient $D_{eff}/\overline{D}_0(\rho _0)$ and mean drift ${V_{eff}} /v_{wall}$ with increasing particle density $\rho _0$ in figures 3(a,b). We find striking variations with increasing density. At low $Pe$, systems of interacting particles differ very little from systems with no interactions: we still observe entropic slow-down. Interestingly, at intermediate $Pe$, the enhancement of the diffusion coefficient increases with particle density $\rho _0$. The turnaround point at a critical value of $Pe$ increases also with increasing particle density. Similarly, mean particle drift is significant at larger values of $Pe$, the more so with increasing particle density. Increasing particle density thus does appear to bridge the gap between the ideal gas case and the incompressible fluid.

Figure 3. Transport properties of a fluid of interacting particles: (a) effective diffusion and (b) mean drift, for a compressible fluid of soft-core interacting Brownian particles with varying density $\rho _0$. Numerical parameters used here are similar to those in figure 2, with additionally $L=200 d_c$, $d_c = 2^{1/6}\ell _0 \simeq 1.12\ell _0$ and $\alpha =1 k_B T/\ell _0^2$. Error bars correspond to one standard deviation over 10 independent runs. Error bars for lower-density values are larger due to a smaller number of tracked particles (see table 1). Theory curves for the ideal gas and incompressible fluid are the same as for figures 2(a,b). Theory curves for the fluid of interacting particles correspond to (4.19) and (4.20).

3.3. Increasing interactions impacts the density profile

To understand the behaviour of $D_{eff}$ and ${V_{eff}}$, we need to understand first how particles rearrange due to inter-particle forces. Typically, a tracer particle will be forced down density gradients because of repulsive pairwise interactions. We therefore measure and analyse the particle density distributions in the channel for various average densities (see figures 4a,b). First, we find that the particle distribution is uniform in the vertical direction, as expected within the lubrication approximation: particles diffuse sufficiently fast on the vertical axis compared to all other relevant time scales. Second, compared to the ideal gas case at the same wall speed (figure 2(c) versus figure 4(a)), the particle distribution is also quite homogeneous on the horizontal direction, much like in the incompressible case (figure 2d). However, at higher wall velocities (figure 4b), particles accumulate again at the bottleneck region. Such accumulation in narrow channels is also observed in simulations of driven, interacting tracers in corrugated channels (Suárez et al. Reference Suárez, Hoyuelos and Mártin2015; Suárez, Hoyuelos & Mártin Reference Suárez, Hoyuelos and Mártin2016), which share a similar geometry in the reference frame where the wall is static.

Figure 4. Particle distribution within the channel for varying particle densities. (a,b) Particle distribution in a fluid of interacting particles at high density $\rho _0=10 \ell _0^{-2}$: (a) for intermediate wall velocity $v_{wall}=0.1 \ell _0/\tau _0$, and (b) for high wall velocity $v_{wall}=0.5 \ell _0/\tau _0$. The colour scale for the density profiles is shared between (a) and (b), and yellow (purple) regions indicate regions of high (low) density. (c,d) Marginal density $p(x)=\int \rho (x,y) \,\text {d} y$ of particles in the channel, calculated by integrating the two-dimensional density profile in (a) over vertical slabs as with the dashed grey box: (c) $p(x)$ for different values of $v_{wall}$, and (d) for different values of particle density $\rho _0$. Numerical parameters used here are similar to those in figure 3. Theory curves correspond to (4.10).

To further quantify the particle distributions, we study the marginal distribution profiles $p (x)=\int _{-\infty }^{\infty }\rho (x,y) \,\text {d} y$ along the channel, obtained in the simulations. As expected, $p(x)$ presents a peak, which indicates particles accumulating near the bottleneck. We find again that the profile $p(x)$ is more peaked with increasing wall velocity (figure 4c). It flattens out with increasing particle density $\rho _0$ (see figure 4d), indicating that density profiles indeed converge to the incompressible fluid case.

4. Analytic theory with pairwise interactions

4.1. Model for particle density profiles along the channel

To understand how the average density profile depends on interaction parameters, we expand our analytical theory. Our starting point is the equation governing the evolution of the density of particles $\rho (x,y,t)$, which can be written formally as (Dean Reference Dean1996; Kawasaki Reference Kawasaki1998)

(4.1)\begin{equation} \frac{\partial{\rho}}{\partial t} = D_0\,\boldsymbol{ \nabla}^2 \rho + \frac{D_0}{k_B T}\,\boldsymbol{ \nabla} \boldsymbol{\cdot} ( \rho\,\boldsymbol{ \nabla} ( \mathcal{V}_{int} * \rho ) ) + \boldsymbol{\nabla} \boldsymbol{\cdot} ( \sqrt{2 D_0 \rho}\,\boldsymbol{ \xi} ),\end{equation}

where $\boldsymbol { \xi }$ is a two-component Gaussian white noise field, such that $\langle \xi _\mu (x,y,t) \rangle = 0$ and $\langle \xi _\mu (x,y,t)\,\xi _\nu (x',y',t') \rangle = \delta _{\mu \nu }\, \delta (x-x')\,\delta (y-y')\,\delta (t-t')$. Here, the convolution product is $(\mathcal {V}_{int} * \rho )(x,y) = \iint \text {d} x' \,\text {d} y'\,\mathcal {V}_{int}(\sqrt {x'^2 + y'^2})\,\rho (x- x', y-y')$, where $\mathcal {V}_{int}$ is any pairwise interaction potential. The term $D_0$ represents the microscopic diffusion constant of the individual Brownian particles. In an infinite domain and in the absence of interactions, the long-time diffusion constant of a tracer will be identical to $D_0$. The effects of interactions and geometry in the problem at hand here both play a role in modifying the late-time diffusion constant with respect to the microscopic one.

To make progress, we first decompose the density into average and fluctuating components:

(4.2)\begin{equation} \rho(x,y,t) = \bar{\rho}(x,y,t)+ \psi(x,y,t), \end{equation}

where $\psi (x,y,t)$ is the noise field such that $\langle \psi (x,y,t) \rangle = 0$. The noise perturbation $\psi (x,y,t)$ can be inferred, as in Démery et al. (Reference Démery, Bénichou and Jacquin2014), by assuming a large, flat, fixed environment and in the case where the background concentration $\bar {\rho }(x,y,t) = \rho _0$ is uniform. The short distance fluctuations in the density field can then be absorbed into a renormalization of the diffusion constant, which now depends on the average density. We will denote by $\overline{D}_0(\rho _0)$ the diffusion coefficient renormalised by the interactions. Here, we will assume that the density fluctuations are short enough in time and space, such that $\bar {\rho }(x,y,t)$ may be treated as a constant, and such that the same renormalisation applies here quickly enough to yield a local diffusion constant depending on the local average density $\overline{D}_0(\bar {\rho })$. This assumption is basically the same as that used in the macroscopic fluctuation theory by Bertini et al. (Reference Bhattacharjee and Datta2015), where it is argued that a coarse-grained hydrodynamic description of a system simply leads to an equation of the type (4.1), but with a local diffusion constant and mobility (that should be coupled to the external forces) that depends on the local density field. Here, we assume that we can write a coarse-grained equation for the evolution of the average density $\bar {\rho }$ as

(4.3)\begin{equation} \frac{\partial{\bar{\rho}}}{\partial t} =- \boldsymbol{ \nabla} \boldsymbol{\cdot}\left( - \overline{D}_0(\bar{\rho})\,\boldsymbol{ \nabla} \bar{\rho} - \frac{\overline{D}_0(\bar{\rho})}{k_B T}\,\bar{\rho}\,\boldsymbol{ \nabla} ( \mathcal{V}_{int} * \bar{\rho} ) \right). \end{equation}

As the field $\bar {\rho }(x,y,t)$ is smooth and the interaction is short-range compared to all the other length scales in the system, we can make the local approximation

(4.4)\begin{equation} \mathcal{V}_{int}(x,y) = \frac{E_0}{\rho_0}\,\delta(x)\,\delta(y), \end{equation}

such that the Fokker–Planck equation for the particle density simplifies to

(4.5)\begin{equation} \frac{\partial{\bar{\rho}}}{\partial t} =- \boldsymbol{ \nabla}\boldsymbol{\cdot} \left( - \overline{D}_0(\bar{\rho})\,\boldsymbol{ \nabla} \bar{\rho} - \overline{D}_0(\bar{\rho})\,\frac{E_0}{k_B T}\,\frac{\bar{\rho}}{\rho_0}\,\boldsymbol{ \nabla} \bar{\rho} \right), \end{equation}

which is nonlinear in $\bar {\rho }$.

Our goal is now to obtain an expression for the average density $\bar {\rho }(x,y,t)$. We seek a solution beyond the trivial mean field (where $\bar {\rho }(x,y,t) = \rho _0$), which makes our approach similar in essence to other derivations of particles on lattices (Illien et al. Reference Illien, Bénichou, Oshanin, Sarracino and Voituriez2018; Rizkallah et al. Reference Rizkallah, Sarracino, Bénichou and Illien2022). Notice that for the height fluctuation profile considered here, the density field is stationary in the frame of reference where the wall is static; i.e. making the change of variables $x' = x - v_{wall}\,t$ and dropping the prime signs yields

(4.6)\begin{equation} 0 = v_{wall}\,\frac{\partial{\bar{\rho}(x,y)}}{\partial x} - \boldsymbol{ \nabla} \boldsymbol{\cdot}\left( - \overline{D}_0(\bar{\rho})\,\boldsymbol{ \nabla} \bar{\rho}(x,y) - \overline{D}_0(\bar{\rho})\,\frac{E_0}{k_B T}\,\frac{\bar{\rho}}{\rho_0}\,\boldsymbol{ \nabla} \bar{\rho} \right). \end{equation}

Integrating vertically, and using exactly the same arguments as in § 2.2.2 for the boundary conditions, shows that

(4.7)\begin{align} J = \int_0^{h(x,t)} \text{d} y \left( - v_{wall}\,\bar{\rho}(x,y) - \overline{D}_0(\bar{\rho})\,\frac{\partial{\bar{\rho}(x,y)}}{\partial x} - \overline{D}_0(\bar{\rho})\,\frac{E_0}{k_B T}\,\frac{\bar{\rho} (x,y)}{\rho_0}\,\frac{\partial{\bar{\rho}(x,y)}}{\partial x} \right) , \end{align}

where $J$, the vertically integrated longitudinal flux, is a constant to be determined.

To solve (4.7), we first assume that we can make the Fick–Jacobs approximation $\bar {\rho }(x,y) \simeq \bar {\rho }(x)$, where the density profile depends only on the longitudinal coordinate, which is correct to first order in $H/L$ and can be demonstrated using lubrication theory. Second, we assume that the average density profile has only small variations around its mean value $\rho _0$, and that the variations are of order $\varepsilon = h_0/H$. We therefore seek a perturbative solution as $\bar {\rho }(x,y) \simeq \rho _0 + \varepsilon \,\rho _1(x) + \cdots$, where $\rho _0 = N/(L H)$ is the mean particle density in the simulation. We can then expand (4.7) in powers of $\varepsilon$. At the lowest order, we find the value of the vertically integrated drift $J = - v_{wall}\,\rho _0 H$. At the next order, we obtain a closed set of equations for the perturbation $\rho _1$,

(4.8)\begin{equation} \left.\begin{gathered} - H\,v_{wall}\,\rho_0 \varepsilon \cos(k_0 x) = H\left[- \overline{D}_0(\rho_0)\left( 1 + \frac{E_0}{k_B T}\right)\frac{\partial{\rho_1(x)}}{\partial x} - v_{wall}\,\rho_1(x) \right] , \\ \int_0^L \rho_1 (x)\,H \,\text{d}x = 0, \\ \rho_1 (x + L) = \rho_1 (x), \end{gathered}\right\} \end{equation}

with relevant periodic boundary conditions, and vanishing integral of $\rho _1$ since the average density should be conserved. We abbreviate $\overline{D}_0 = \overline{D}_0(\rho _0)$ in the following.

Interestingly, we find in (4.8) that the density profile relaxes with an effective diffusion coefficient

(4.9)\begin{equation} D^\star \equiv \overline{D}_0 \left( 1 + \frac{E_0}{k_B T}\right)= \overline{D}_0\,\frac{\chi_T^{id}}{\chi_T^{{\star}}} \end{equation}

that is the collective diffusion coefficient, as it characterises how the fluid of interacting particles relaxes, not how a single particle diffuses. Equation (4.9) is also expected for hard spheres (Hess & Klein Reference Hess and Klein1983; Lahtinen et al. Reference Lahtinen, Hjelt, Ala-Nissila and Chvoj2001). If interactions are repulsive (i.e. $E_0>0$), then the fluid becomes more incompressible with $\chi _T^\star <\chi _T^{id}$, and density inhomogeneities relax faster than they would in the ideal gas, since $D^\star >\overline{D}_0$. Conversely, attractive interactions (i.e. $E_0<0$) increase compressibility and stabilise density inhomogeneities. This property is absolutely essential to understand the long-time transport properties in a fluctuating medium.

Solving (4.8), we find

(4.10)\begin{equation} \bar{\rho}(x) = \rho_0 + \varepsilon\,\rho_1(x) = \rho_0\left( 1 - \frac{h_0}{H}\,\frac{Pe^{{\star}}}{1 + (Pe^{{\star}})^2}\left[ \sin(k_0x) + Pe^{{\star}} \cos(k_0x) \right] \right), \end{equation}

and a new and natural Péclet number characterising the fluctuations emerges,

(4.11)\begin{equation} Pe^{{\star}} = \frac{v_{wall}}{D^\star k_0} = \frac{\omega_0}{D^\star k_0^2}, \end{equation}

and takes into account the enhanced collective diffusion $D^\star$ due to repulsive interactions. Notice that the solution can be rewritten further as

(4.12)\begin{equation} \bar{\rho}(x) = \rho_0 + \varepsilon\,\rho_1(x) = \rho_0\left( 1 + \frac{h_0}{H}\,\frac{Pe^{{\star}}}{\sqrt{1 + (Pe^{{\star}})^2}} \cos(k_0x + \varphi) \right), \end{equation}

where $\varphi = {\rm \pi}- \arctan ( 1/Pe^{\star } )$. This shows that the perturbation in the density of particles due to the fluctuating interface propagates with a phase $\varphi$ that is characterised by how fast particles relax as a group.

We plot the resulting longitudinal probability density profiles $p(x) = [\rho _0 + \varepsilon \, \rho _1(x)]/L\rho _0$, where $\rho _1(x)$ is given by (4.10), in figures 4(c,d), and find excellent agreement with the numerical results. More specifically, particles accumulate at the constriction at high forcings $v_{wall}$ (figure 4c, light purple) or at low particle densities (figure 4d, yellow and orange). If collective effects are weak, then the time scale for density fluctuations to diffuse across a bulge $1/D^\star k_0^2$ is larger than the time that it takes a bulge to move, $1/v_{wall}\,k_0 = 1/\omega _0$, or equivalently, $Pe^{\star } \gg 1$. Clearly, from (4.12), the phase is $\varphi \simeq {\rm \pi}$ and the interface squeezes particles exactly out of phase, i.e. in the constriction. This can also be understood in the conservation of mass (4.8), where at large $v_{wall}$, we simply have $(H - h(x))\,v_{wall}\,\rho _0 = \varepsilon \,\rho _1(x)\,v_{wall}$ such that $\rho _1(x) \sim 1/h(x)$. Similarly as for the isolated particles in § 2, this is a ‘traffic jam’ effect: in the referential where the wall is fixed, particles trying to move with the fast velocity $-v_{wall}$ accumulate where the road is narrow, i.e. at the constriction.

However, there is a trade-off between accumulation due to wall speed, and increased local repulsion of particles in the accumulation region. We observe that at higher particle densities $\rho _0$, the accumulation effect is tempered (figure 4d, orange to purple). When particle repulsion increases (due to either increased particle density $\rho _0$ or interaction strength $\alpha >0$), the local repulsion dominates, resisting compression of the gas in the narrow constriction. This is coherent with the fact that the compressibility $\chi _T^{\star }/\chi _T^{id}$ decreases with particle density $\rho _0$ (or with the strength of the repulsive interactions $\alpha$). As a result, collective effects are strong, $Pe^{\star } \ll 1$, density fluctuations relax faster than the interface's motion, and the density profile flattens out.

We now turn to understand how this local distribution of particles affects the long-time transport properties of a tracer.

4.2. Effective diffusion and drift with interacting particles

4.2.1. Equation for the diffusion of the tracer

We now consider the motion of a tracer in this gas of soft-core interacting Brownian particles, and we infer how the confined fluctuating environment modifies the long-time tracer motion. Based on reasoning similar to that given above, the probability distribution $\rho _{tr}(x,y,t)$ of finding the tracer particle at position $(x,y)$ at time $t$ satisfies

(4.13) \begin{equation} \frac{\partial{\rho_{tr}(x,y,t)}}{\partial t} =- \boldsymbol{ \nabla} \boldsymbol{\cdot} \left(- \overline{D}_0(\bar{\rho})\,\boldsymbol{\nabla} \rho_{tr} - \overline{D}_0(\bar{\rho})\, \frac{E_0}{k_B T}\,\frac{\rho_{tr}}{\rho_0}\,\boldsymbol{ \nabla} \bar{\rho} \right). \end{equation}

Again looking at the longitudinal transport within the Fick–Jacobs framework, the equation on the marginal distribution $p_{tr}(x,t) = \int _0^{h(x,t)} \rho _{tr}(x,y,t) \,\textrm {d} y$ is given by (2.10):

(4.14) \begin{equation} \frac{\partial{p_{tr}(x,t)}}{\partial t} =- \frac{\partial{}}{\partial x} \left(- \overline{D}_0(\bar{\rho})\,\frac{\partial{p_{tr}}}{\partial x} + \overline{D}_0(\bar{\rho})\,p_{tr}\,\frac{\partial{\ln h(x,t)}}{\partial x} - \overline{D}_0(\bar{\rho})\,\frac{E_0}{k_B T}\,\frac{p_{tr}}{\rho_0}\,\frac{\partial{\bar{\rho}}}{\partial x} \right). \end{equation}

To understand how the tracer's motion is modified by interactions with other particles and by the fluctuating boundary, we expand (4.14) in powers of $\varepsilon =h_0/H$, and obtain the particle's long-time diffusion coefficient and drift, at order $\varepsilon ^2$, as was done in § 2.3. Since the average particle density depends on space, the local diffusion coefficient of the tracer also depends on space as $\overline{D}_0(\bar {\rho }(x,t)) = \overline{D}_0(\rho _0) + \varepsilon \,\overline{D}_0'(\rho _0)\,({\rho _1(x,t)}/{\rho _0}) + O(\varepsilon ^2)$. Here, since the prescribed wall fluctuations $h(x,t)$ are periodic in space, it is possible to obtain implicit expressions for the long-time effective diffusion coefficient $D_{eff}$ and drift $V_{eff}$ using the approach reported in Reimann et al. (Reference Reimann, Van den Broeck, Linke, Hänggi, Rubi and Pérez-Madrid2001) and in Guérin & Dean (Reference Guérin and Dean2015), and resolve their explicit expressions after cumbersome expansions in the small parameter $\varepsilon$ (see § 3 of the supplementary material). In the cases explored here, the local change in diffusion coefficient is small, $\overline{D}_0'(\rho _0)/\overline{D}_0(\rho _0)\,\rho _0 \ll 1$, and we find that its impact on the tracer dynamics can be neglected. We can therefore use the explicit perturbation theory results in § 2.3, assuming $\overline{D}_0(\bar {\rho }(x,t)) \simeq \overline{D}_0(\rho _0) \equiv \overline{D}_0$. As expected, the perturbation theory and the periodic framework approach of Reimann et al. (Reference Reimann, Van den Broeck, Linke, Hänggi, Rubi and Pérez-Madrid2001) provide exactly the same results, to leading order in $\varepsilon$. For the sake of completeness, we will nonetheless derive, in a future work, a general perturbation theory with explicit results for Fokker–Planck equations with diffusion and drift with arbitrary dependence on space and time in any dimension.

4.2.2. Local drift of a tracer in a bath of interacting particles

The equation of motion (4.14) can be simplified to

(4.15)\begin{equation} \frac{\partial{p_{tr}(x,t)}}{\partial t} =- \frac{\partial{}}{\partial x} \left(- \overline{D}_0\,\frac{\partial{p_{tr}(x,t)}}{\partial x} + v(x,t)\,p_{tr}(x,t) \right), \end{equation}

where the local longitudinal drift is

(4.16)\begin{equation} v(x,t) = \overline{D}_0\,\frac{1}{H}\,\frac{\partial{h(x,t)}}{\partial x} - \overline{D}_0\,\frac{E_0}{k_B T}\,\frac{h_0}{H}\,\frac{1}{\rho_0}\, \frac{\partial{\rho_1(x,t)}}{\partial x} \end{equation}

at lowest order in $\varepsilon = h_0/H$. We can verify that terms of $O(\varepsilon ^2)$ in the local drift vanish when averaged over one period and hence do not contribute to the renormalization of the long-time transport properties. The last term of (4.16) quantifies the effect of interactions on local tracer motion, and was not present in the ideal gas case. It indicates that particles drift away from accumulation regions because of repulsive interactions. Finally, the magnitude of the effect is proportional to the interaction strength $E_0/k_B T = (\chi _T^{id} - \chi _T^{\star })/\chi _T^{\star }$, or to the relative compressibility of the fluid.

We compute in simulations the mean local particle velocity in the interacting gas of particles; see figure 5. With increasing incompressibility (increasing $\rho _0$ in figures 5ac), a local drift emerges that carries particles away from the accumulation region located at the channel constriction. As expected with increasing interactions, the velocity field approaches that of the incompressible fluid (reported for comparison in figure 5d, already shown in figure 2d).

Figure 5. Local drift and particle density profiles within the channel for varying mean particle densities. (ac) Simulation results for increasing mean density for simulations with interacting particles. (d) Simulation results for the incompressible fluid case. Note that (d) corresponds to figure 2(d) and is shown as a reminder. The colour scale for the density profiles is shared between all plots. Yellow (purple) regions indicate regions of high (low) density. The blue arrows correspond to the local velocity of particles (averaged over time), and the length of the arrows is proportional to the magnitude of the velocity. The scale of the arrows is arbitrary but is the same across all plots, hence small arrows in (a) indeed indicate very weak velocity fields compared to (c) or (d). The interaction strength for (ac) is $\alpha = 1 k_B T/\ell _0^2$, and other numerical parameters are similar to those in figure 3. For all plots, $v_{wall} = 0.5 \ell_0 /\tau_0$.

These results are confirmed by the analytical prediction, using the expression for $\rho _1$ in (4.10). Recalling that the collective Péclet number is $Pe^{\star } = \omega _0^2 / D^\star k_0$, and denoting $\bar {Pe}=\omega _0^2 / \overline{D}_0 k_0$ as the Péclet number of a tracer in medium of density $\rho _0$, we obtain in Fourier space

(4.17)\begin{equation} \tilde{v}(k,\omega) = \frac{\tilde{h}}{H}\,\frac{{\rm i} \overline{D}_0 k (1 + Pe^{{\star}} \bar{Pe}) + (Pe^{{\star}}/\bar{Pe} - 1)\omega/k }{{1 + (Pe^{{\star}})^2}}. \end{equation}

Equation (4.17) interpolates perfectly between the ideal gas and the incompressible fluid. When interactions vanish, $D^\star =\overline{D}_0=D_0$ and $Pe^{\star }=\bar {Pe}=Pe$, which yields the ideal gas expression for the local velocity (2.16), where $\tilde {v}(k,\omega ) = \textrm {i} D_0 k \tilde {h}/H$. Reciprocally, the incompressible fluid limit is $D^\star \to \infty$, which implies $Pe^{\star } \rightarrow 0$, hence ${\tilde {v}(k,\omega ) = (\textrm {i} \overline{D}_0 k - \omega /k)\tilde {h}/H}$, as found already in (2.22) but with $\overline{D}_0$ instead of $D_0$ since the interactions have renormalised the bare diffusion coefficient of the tracer.

4.2.3. Effective long-time diffusion and mean drift

We use the perturbation results in (2.14) and (2.15) to obtain the long-time diffusion coefficient and drift. From (4.17), the spectrum of the local velocity fluctuations is

(4.18)\begin{equation} S(k,\omega) = \frac{\overline{D}_0^2 k^2}{(2{\rm \pi})^2}\,\frac{1 + \bar{Pe}^2}{1 + (Pe^{{\star}})^2}\,\frac{h_0^2}{4 H^2} \left( \delta(k+k_0)\,\delta(\omega - \omega_0) + \delta(k-k_0)\,\delta(\omega + \omega_0) \right), \end{equation}

and the long-time transport properties are

(4.19)\begin{gather} \frac{D_{eff}^{\textit{int. gas}}}{ \overline{D}_0} = 1 - \frac{h_0^2}{2H^2}\,\frac{1 - 3 \bar{Pe}^2}{1 + \bar{Pe}^2}\,\frac{1}{1 + (Pe^{{\star}})^2}, \end{gather}
(4.20)\begin{gather}\frac{V_{eff}^{\textit{int. gas}}}{v_{wall}} = \frac{h_0^2}{2H^2}\,\frac{1}{ 1 + (Pe^{{\star}})^2}. \end{gather}

Equations (4.19) and (4.20) form the main analytic result of the paper and give the transport properties of tracers in a fluid with arbitrary compressibility. We recall that the modified Péclet number is connected to the compressibility $Pe^{\star } = \bar {Pe} \chi _T^{\star }/\chi _T^{id}$, where $\bar {Pe}$ is the Péclet number characterising the fluctuations. For weak interactions, where $\rho _0, \alpha \rightarrow 0$, the compressibility is that of the ideal gas $\chi _T^{\star } \rightarrow \chi _T^{id}$, and $Pe^{\star } \rightarrow Pe$. In this limit, we recover the ideal gas results (2.18) and (2.19). Reciprocally, for strong interactions, $\rho _0, \alpha \rightarrow \infty$, $\chi _T^{\star } \rightarrow 0$, hence $Pe^{\star } \rightarrow 0$ and we recover the incompressible fluid results of (2.24) and (2.25). Note that in the limit $\rho _0\to \infty$, $\overline{D}_0(\rho _0)$ is still finite; see § A.4. Equations (4.19) and (4.20) therefore interpolate between the ideal gas and the incompressible fluid cases.

Despite the mean-field-like assumptions made, our analytic theory summarised in (4.19) and (4.20) perfectly captures the simulation results (see the solid lines in figures 3 and 10), and confirms transport mechanisms in this complex environment. With increasing particle density, particle collisions push particles away from the accumulation region, further into the bulges. Particle–wall collisions then push and disperse particles. In this complex environment, and in contrast with the paradigm where crowded environments slow down diffusion (Lekkerkerker & Dhont Reference Lekkerkerker and Dhont1984; Lowen & Szamel Reference Lowen and Szamel1993; Dean & Lefèvre Reference Dean and Lefèvre2004), here particle collisions or interactions favour mixing.

5. Discussion and conclusion

In this work, we have explored the impact of crowding on the long-time transport properties of particles in fluctuating channels. Our simulation results show a broad range of behaviours that are well captured by our explicit analytic theory based on a perturbation expansion in the wall fluctuation amplitude $h_0/H$. The results are best described in terms of a Péclet number characterising the fluctuations, $Pe= \omega _0/D_0 k_0^2 = 2{\rm \pi} v_{wall}/D_0 L$. At low $Pe \ll 1$, corresponding to fixed interfaces, all fluids behave similarly: particles are slowed down by constrictions, and the effective diffusion is decreased. At intermediate $Pe \simeq 1$, particle–wall collisions stir particles, and the effective diffusion is increased. This effect persists until the wall moves so fast that particles no longer have time to reach the moving bulges and accumulate in the constrictions. In this final regime, $Pe \gg 1$, the effective diffusion is unchanged. The accumulation regime arises for higher and higher $Pe$ values for increasing particle–particle interactions, i.e. for increasing incompressibility that resists accumulation.

5.1. Collisions enhance diffusion

One of the main findings of our work is that here, both numerically and analytically, we have demonstrated that increasing repulsive interactions or collisions between particles can enhance the late-time diffusion coefficient and the mean drift characterising the dispersion of a tracer particle in fluctuating channels. This is in contrast with the intuition that collisions in equilibrium reduce the diffusion coefficient (Lekkerkerker & Dhont Reference Lekkerkerker and Dhont1984; Lowen & Szamel Reference Lowen and Szamel1993). The mechanism of diffusion enhancement is in fact rather simple: collisions or repulsive interaction help to push particles closer to the walls. Eventually, wall–particle collisions help mixing. Since this mechanism is rather straightforward, we expect it to be broadly applicable – for example, beyond lubrication approximation or for hard-core repulsive interactions. Such effects could be explored within our framework or alternatively using dynamic density functional theory (Marconi & Tarazona Reference Marconi and Tarazona1999). In more detailed physical settings, such as with hard-core interactions, other effects would likely also come into play; for example, the accessible volume in the channel is smaller for larger particles (Riefler et al. Reference Riefler, Schmid, Burada and Hänggi2010; Suárez et al. Reference Suárez, Hoyuelos and Mártin2015), and hydrodynamic effects become important (Yang et al. Reference Yang, Liu, Li, Marchesoni, Hänggi and Zhang2017).

We now put our results in a broader context. In the Introduction, we recalled that diffusion of a driven, out-of-equilibrium tracer in a bath of interacting particles is enhanced by repulsive interactions (Bénichou et al. Reference Bénichou, Illien, Oshanin and Voituriez2013, Reference Bertini, De Sole, Gabrielli, Jona-Lasinio and Landim2018; Démery et al. Reference Démery, Bénichou and Jacquin2014; Illien et al. Reference Illien, Bénichou, Oshanin, Sarracino and Voituriez2018). In a confined channel, thermal fluctuations of the wall could possibly enhance the diffusion coefficient of particles, as we have seen in the limiting case of an incompressible fluid (Marbach et al. Reference Marbach, Dean and Bocquet2018). Another recent work finds that diffusion of odd-diffusing interacting particles is enhanced with increasing densities (Kalz et al. Reference Kalz, Vuijk, Abdoli, Sommer, Löwen and Sharma2022). While the physical set-up in Kalz et al. (Reference Kalz, Vuijk, Abdoli, Sommer, Löwen and Sharma2022) is very different from ours, the mathematical similarities that lead to diffusion enhancement are striking (for example, comparing their Eq. (9) with our (2.24)), and one might speculate that there exists a universal framework to understand these effects under the same light.

5.2. Open questions on fluctuating channels

Beyond the question asked here, namely of understanding how crowding affects transport in fluctuating channels, there are many open fundamental questions. For example, boundaries are not necessarily repulsive and smooth. Surface rugosity would lead to Taylor dispersion in an incompressible fluid (Marbach & Alim Reference Marbach and Alim2019; Kalinay Reference Kalinay2020), but how would surface rugosity of the wall potential induce Taylor dispersion in the fluid of interacting particles? Attraction at the boundaries also leads to surprising speed-up of diffusion in corrugated static channels (Alexandre et al. Reference Alexandre, Mangeat, Guérin and Dean2022). Is this speed-up enhanced further by fluctuations? Considering more complex fluids would likely modify the distribution of particles along the channel, and hence the transport properties. For example, normal stresses close to the boundaries in odd viscous fluids bring particles close to the walls (Hargus et al. Reference Hargus, Klymko, Epstein and Mandadapu2020; Fruchart, Scheibner & Vitelli Reference Fruchart, Scheibner and Vitelli2023). Could this phenomenon enhance dispersion? Down the scales, molecular (Yoshida et al. Reference Yoshida, Kaiser, Rotenberg and Bocquet2018) or quantum (Kavokine, Bocquet & Bocquet Reference Kavokine, Bocquet and Bocquet2022; Coquinot, Bocquet & Kavokine Reference Coquinot, Bocquet and Kavokine2023) effects enhance the mobility of individual molecules; how would these effects combine with mechanical fluctuations? With the advent of highly sensitive techniques to probe the motion of particle near surfaces in soft and increasingly complex environments (Zhang et al. Reference Zhang, Bertin, Arshad, Raphael, Salez and Maali2020; Sarfati et al. Reference Sarfati, Calderon and Schwartz2021; Vilquin et al. Reference Vilquin, Bertin, Raphaël, Dean, Salez and Mcgraw2022), one might expect to answer some of these questions in the light of further experimental results.

Supplementary material

Supplementary material is available at


The authors are grateful for fruitful discussions with L. Bocquet, J.-P. Hansen, P. Levitz, M.-T. Hoang Ngoc, G. Pireddu, B. Rotenberg, B. Sprinkle and A. Thorneywork. S.M., R.Z. and Y.W. thank the Applied Math Summer Undergraduate Research Experience Program of the Courant Institute for putting them in contact. S.M. would like to thank the Institut d’Études Scientifiques de Cargése for hosting the Transport in Narrow Channels workshop that led to inspiring discussions for this work. R.Z. would like to thank the Laboratoire MSC Paris and the Center for Data Science ENS Paris for hospitality.


This work was supported in part by the MRSEC Program of the National Science Foundation under award number DMR-1420073. S.M. was also supported by the European Union's Horizon 2020 research and innovation program under the Marie Skłodowska-Curie grant agreement 839225 MolecularControl. R.Z. was also supported by grant no. NSF DMR-1710163.

Declaration of interests

The authors report no conflict of interest.

Data availability statement

The data that support the findings of this study are available upon reasonable request to the author.

Appendix A. Simulation details

A.1. Simulation algorithm

All simulations are performed using a forward Euler stochastic scheme to discretise the overdamped Langevin dynamics of the particles. For a particle $i$ at position ${\boldsymbol {X}_i(t) = (x_i(t),y_i(t))}$, the following position at time $t + \Delta t$ is computed as

(A1)\begin{equation} \boldsymbol{X}_i(t+\Delta t) = \boldsymbol{X}_i(t) + \Delta t\,\boldsymbol{U}_i + \Delta t\,\frac{D_0}{k_B T}\,(\boldsymbol{F}_i(t))+\sqrt{2D_0\,\Delta t}\,\boldsymbol{G}_i(t), \end{equation}

where $\boldsymbol {U}_i = \boldsymbol {U} (x_i(t),y_i(t))$ is the background flow field (which is non-zero only in the incompressible case), $\boldsymbol {F}_i = \boldsymbol {F}_{wall} + \sum _{j\neq i} \boldsymbol {F}_{{int}, ij}$ is the sum of the forces exerted by the channel walls, $\boldsymbol {F}_{wall}$, on the particle and by the neighbouring particles $\boldsymbol {F}_{{int}, ij}$ when interactions are present, and $\boldsymbol {G}_i(t)$ is a vector of two independent random numbers drawn from a Gaussian distribution of mean 0 and variance 1. Note that our simulations have been tested with a time step twice as small and yielded no significant difference.

Particles are confined in the channel by means of a potential that exerts a force on the particles only if they reach the boundaries. More precisely, the force exerted by the wall on a particle with coordinate $(x,y)$ is given by $\boldsymbol {F}_{wall}=-\boldsymbol { \nabla } \mathcal {V}_{wall}(x,y)$, with

(A2)\begin{equation} \mathcal{V}_{wall}(x,y) =\begin{cases} \lambda (y-Y_{upper}(x,t))^4 & \text{if } y>Y_{upper}(x,t),\\ 0 & \text{if } Y_{lower}(x,t)< y< Y_{upper}(x,t),\\ \lambda (y-Y_{lower}(x,t))^4 & \text{if } y< Y_{lower}(x,t), \end{cases} \end{equation}

with $\lambda >0$ a stiffness coefficient (with dimensions of energy over length to the power 4), and where the boundary equations are given by

(A3)\begin{gather} Y_{upper}(x,t) = H + h_0 \cos\left( \frac{2{\rm \pi}}{L}\,(x-v_{wall}\,t)\right), \end{gather}
(A4)\begin{gather}Y_{lower}(x,t) =0, \end{gather}

where $H$ and $h_0$ represent the average height and the variation amplitude of the channel height, respectively. Using a soft confining potential to model the wall is convenient for simulation purposes, as it avoids dealing with reflecting Brownian walks, which carries some challenges (Scala, Voigtmann & De Michele Reference Scala, Voigtmann and De Michele2007). It also allows one to keep a rather large integration time step. Finally, this choice of boundary conditions bears a physical meaning since it considers that the wall is likely made of a slightly penetrable material that would interact with the particles/colloids, and is, therefore, more general than a hard-reflecting boundary condition. Such boundary models have been used extensively in theoretical active matter systems (Solon et al. Reference Solon, Fily, Baskaran, Cates, Kafri, Kardar and Tailleur2015; Zakine et al. Reference Zakine, Zhao, Knežević, Daerr, Kafri, Tailleur and van Wijland2020; Ben Dor et al. Reference Ben Dor, Ro, Kafri, Kardar and Tailleur2022), and as our theory relies only on the presence of a boundary, the results are largely unaffected by a change of potential, as long as the boundary layer of particles at the wall and subjected to the potential repulsion is small compared to all other relevant length scales in the system. The fact that our numerical results (and especially the two-dimensional density profiles) are consistently in agreement with hard-reflecting boundary conditions taken in the analytical computation shows that this choice of boundary conditions is fairly equivalent to hard-reflecting boundaries.

In § 3, we have performed simulations with pairwise interacting particles. The force on particle $i$ exerted by particle $j$ is $\boldsymbol {F}_{{int}, ij} = - \boldsymbol { \nabla } \mathcal {V}_{int}(r_{ij})$, where $r_{ij}$ is the distance between $i$ and $j$. We chose simple repulsive interactions described by the soft potential (3.1), which we recall here as $\mathcal {V}_{int}(r) = \alpha (r - d_c)^2\,\varTheta (r < d_c)$, where $\alpha$ is the interaction strength in units of a spring constant, and $d_c$ is the typical particle diameter (see figure 6). For this quadratic interaction $\mathcal {V}_{int}$, the expression of its vertically integrated version is given by

(A5)\begin{align} \mathcal{U}_{int}(x) &= \int_{-\infty}^{+\infty} \mathcal{V}_{int}(x,y) \,\text{d} y \end{align}
(A6)\begin{align} &= 2 \alpha \left(\frac{1}{3}\,(d_c^2+2 x^2) \sqrt{d_c^2-x^2} +d_c x^2 \log \left(\frac{x}{\sqrt{d_c^2-x^2}+d_c}\right)\right) \varTheta(d_c > |x|). \end{align}

This expression is used to compute the density profiles $\rho (x)$, but we always use the radial potential $\mathcal {V}_{int}$ in the Monte Carlo simulations.

Figure 6. Radial potential $\mathcal {V}_{int}$ (blue) and its one-dimensional smoothed out version $\mathcal {U}_{int}$ (red). Parameters are $d_c=1 \ell _0$ and $\alpha =2 k_B T / \ell _0^2$.

The system size along the $x$ direction is $L$, and we take periodic boundary conditions in the $x$ direction. We keep track of the full dynamics of the $N$ particles with their absolute position $\boldsymbol {X}_{true}$ to compute their mean-square displacement and mean drift, but interactions are computed using the folded positions in the periodic domain $[0,L)\times (-\infty, \infty )$. The mean density of particles is $\rho _0 = N/(L H)$. To access higher particle densities, we therefore change the total number of particles $N$ in the simulation (see table 1).

Table 1. Typical simulation parameters for an interacting particle system.

We initialise systems from a uniform distribution of particles in the bottom part of the channel ($0< y < H- h_0$). We first let the system equilibrate for a time $t_{eq}\sim H/D_0\sim 100 \tau_0$ within a fixed undulated channel ($v_{wall} = 0$). Then the actual simulation starts with a positive value of $v_{wall}$.

A.2. Simulation parameters

The domain characteristics are $L=200\ell _0$, $H=12\ell _0$, $h_0=3\ell _0$ for all the simulations. To simulate non-interacting particles, we take wall stiffness $\lambda =300 k_B T/\ell _0^4$, and the integration time step is $\Delta t=4\times 10^{-3} \tau _0$. For a typical simulation, the total number of iterations is $N_{iter}=5\times 10^6$, and the number of particles is $N=10^5$. To simulate interacting particles, we take $\lambda =10 k_B T/\ell _0^4$, and the integration time step is $\Delta t=4\times 10^{-3} \tau _0$. Typical values of the number of particles $N$ and total number of iterations $N_{iter}$ are shown in table 1.

A.3. Simulation analysis

We perform multiple independent simulations to increase statistical resolution. For each value of $v_{wall}$, we perform 10 independent simulations starting from different initial configurations (and different seeds). Symbols in all graphs represent the mean value of the observable, and error bars correspond to one standard deviation over these 10 independent measurements. The mean drift is calculated simply as $V_{eff} = ({1}/{N})\sum _{i=1}^N \langle {(x_{{true}, i} (t) - x_{true}(0))}/{t} \rangle _t$. The mean-square displacement of particles is calculated as ${MSD}(t) = ({1}/{N})\sum _{i=1}^N \langle (x_{{true}, i} (t + t_0) - x_{true}(t_0))^2 \rangle _{t_0}$, where the average is done over initial times $t_0$. The self-diffusion coefficient is then obtained as a least squares linear fit of the mean-square displacement. The parameters that we choose allow us to neglect the finite size corrections due to periodic boundary conditions. Indeed, such corrections scale as ${\sim d_c/L=0.005}$ (Dünweg & Kremer Reference Dünweg and Kremer1991, Reference Dünweg and Kremer1993) or as $\sim (H/L)^2=0.004$ (Simonnin et al. Reference Simonnin, Noetinger, Nieto-Draghi, Marry and Rotenberg2017), thus are negligible in the numerical measurements performed here.

A.4. Simulation calibration: self-diffusion coefficient $\overline{D}_0(\rho _0)$ of soft-core interacting particles

To calibrate our model, we calculate the long-time self-diffusion coefficient of particles in a fluid of soft-core interacting particles. In figure 7, we report the results of the computation of the diffusion coefficient $\overline{D}_0(\rho _0)$ for fixed flat walls according to different parameters of the interaction: the particle density $\rho _0$, and the interaction strength $\alpha$.

Figure 7. Effective diffusion of a tracer in a homogeneous bath of soft-core interacting particles as a function of (a) the particle density $\rho _0$, and (b) the interaction strength $\alpha$. The tracer is identical to the particles of the bath. Solid lines are obtained from the computation in Démery et al. (Reference Démery, Bénichou and Jacquin2014), and symbols correspond to simulation results. For (a), we fix $\alpha = 1 k_BT/\ell _0^2$ and channel height $H=20\ell _0$, and we simulate for $N=[1,2,5,10,20,10,20] \times 10^3$ particles in a two-dimensional flat channel of length $L=[500,500,500,500,500,100,100] \ell _0$ with periodic boundary conditions. For (b), a couple of values of $\rho _0$ (indicated in the legend) are used. Error bars correspond to one standard deviation over 10 independent runs.

As a self-consistent check, we also compare our simulations to the analytic mean-field theory of Démery et al. (Reference Démery, Bénichou and Jacquin2014), and find good agreement between simulation and theory. We recall briefly the analytic formula here. The correction to the bare diffusion coefficient is given by

(A7)\begin{equation} \frac{\overline{D}_0(\rho_0)}{D_0} \simeq 1-\frac{1}{2 d \rho_0 }\int \frac{\left(\dfrac{\rho_0\,\tilde{\mathcal{V}}_{int}(\boldsymbol{k})}{k_B T}\right)^2}{\left[1+\dfrac{\rho_0\,\tilde{\mathcal{V}}_{int}(\boldsymbol{k})}{k_B T}\right] \left[1+\dfrac{\rho_0\, \tilde{\mathcal{V}}_{int}(\boldsymbol{k})}{2 k_B T} \right]}\,\frac{\text{d}^d\boldsymbol{k}}{(2{\rm \pi})^d}, \end{equation}

where $\rho _0=N/V$ is the particle density, $d$ is the number of spatial dimensions, and $\tilde{\mathcal {V}}_{int}(\boldsymbol {k})$ is the Fourier transform of the interaction pair potential, here

(A8)\begin{equation} \tilde{\mathcal{V}}_{int}(\boldsymbol{k}) = \int \text{d}^d \boldsymbol{r} \,{\rm e}^{-{\rm i} \boldsymbol{k} \boldsymbol{\cdot} \boldsymbol{r}}\,\mathcal{V}_{int}(\boldsymbol{r}). \end{equation}

Note that at large densities $\rho _0$ or interaction strengths, (A7) is no longer valid, and it has been argued that it should be approached by (Dean & Lefèvre Reference Dean and Lefèvre2004)

(A9)\begin{equation} \frac{\overline{D}_0(\rho_0)}{D_0} \simeq \exp \left( -\frac{1}{2 d \rho_0 }\int \frac{\left(\dfrac{\rho_0\, \tilde{\mathcal{V}}_{int}(\boldsymbol{k})}{k_B T}\right)^2}{\left[1+\dfrac{\rho_0\,\tilde{\mathcal{V}}_{int}(\boldsymbol{k})}{k_B T}\right] \left[1+\dfrac{\rho_0\,\tilde{\mathcal{V}}_{int}(\boldsymbol{k})}{2 k_B T} \right]}\,\frac{\text{d}^d\boldsymbol{k}}{(2{\rm \pi})^d} \right). \end{equation}

In our set-up for $d=2$ and $k=|\boldsymbol {k}|$, we have

(A10)\begin{equation} \tilde{\mathcal{V}}_{int}(\boldsymbol{k}) = \frac{2 {\rm \pi}\alpha d_c^2 (-{\rm \pi}\,H_1(d_c k)\,J_0(d_ck)+{\rm \pi}\,H_0(d_c k)\,J_1(d_c k)+2 J_0(d_c k)-4 J_1(d_c k)/(d_c k))}{k^2}, \end{equation}

where $J_\nu (z)$ is the special Bessel function, and $H_\nu (z)$ is the Struve function

(A11)\begin{equation} H_\nu(z) = \left( \frac{z}{2}\right)^{\nu+1}\sum_{n=0}^\infty \frac{(-1)^n (z/2)^{2n}}{\varGamma(n+3/2)\,\varGamma(n+\nu+3/2)}. \end{equation}

We plot (A7) in figure 7, and compare it to our numerical results. Unsurprisingly, the mean-field approximation starts to fail as the interaction strength $\alpha$ increases and prevents particles from overlapping. To present the results of the long-time diffusion coefficients $D_{eff}$ in confined wiggling spaces relative to $D_0$ (figures 3 and 10), we always use the numerically calculated diffusion coefficient $\overline{D}_0(\rho _0)$ in the flat fixed space.

In the main text, we investigate limit behaviours as $\rho _0 \rightarrow \infty$ and $\alpha \rightarrow \infty$. A convenient Gaussian interaction potential can be considered to gain analytical insights into the diffusion. With $\mathcal {V}_{int}(\boldsymbol {r})=\alpha \ \textrm {e}^{-\boldsymbol {r}^2/(2d_c^2)}$, we have $\tilde {\mathcal {V}}_{int}(\boldsymbol {k})=2{\rm \pi} \alpha d_c^2 \ \textrm {e}^{-d_c^2\boldsymbol {k}^2}$in dimension $d=2$, for instance. The integral in (A9) can thus be approximated when integrating only in the range ${\rho _0\,\tilde {\mathcal {V}}_{int}(\boldsymbol {k})}/{(k_B T)}\gg 1$. The interaction-and-density-dependent cutoff is given by $\varLambda (\rho _0)\sim ({1}/{d_c})\sqrt {\ln (\alpha d_c^2 \rho _0/k_B T)}$, and the diffusion is thus given by

(A12)\begin{equation} \frac{\overline{D}_0(\rho_0)}{D_0} \simeq \exp \left( -\frac{(2\varLambda)^d}{d \rho_0 (2{\rm \pi})^d } \right), \end{equation}

which yields, first, $\overline{D}_0(\rho _0)\to D_0$ for $\rho _0\to \infty$ and $\alpha$ fixed (as expected since the potential landscape becomes flat), and second, $\overline{D}_0(\rho _0)\to 0$ for $\alpha \to \infty$ and $\rho _0$ fixed (as expected for jamming).

A.5. Simulation calibration: confinement with the soft wall potential

For each simulation type, we check the penetration depth of the particles in the confining wall. With increasing repulsive interactions (with increased $\rho _0$ or $\alpha$), and since the wall is ‘soft’, particles may be squeezed into the confining soft wall. We then estimate the penetration depth of each fluid within the confining wall by looking at the probability density at the centre of the channel where the constriction is (in the frame of reference where the wall is fixed); see figure 8. We find that indeed, with increasing interactions the penetration depth $\delta h = (0.2- 1) \ell _0$ increases. We record the penetration depth $\delta h$ from figure 8 for each set of numerical parameters, corresponding to the depth for which the probability density is half of its bulk value (dashed black horizontal line). We then use $H = H + 2 \delta h$ (since there is an upper wall and a bottom wall) in all analytical formulas.

Figure 8. Vertical density profile near the repulsive confining wall of particle systems for (a) various particle densities and (b) various interaction strengths. Numerical systems correspond to those detailed in figures 3 and 10, and the confining wall is set in both cases as $y = h(x = L/2) = H - h_0 = (12 - 3)\ell_0 = 9\ell _0$. Hence particle systems extend up to $\delta h \sim \ell _0$ into the confining wall.

Appendix B. Additional data for the ideal gas

To test the validity of the analytic derivation in § 2 for the transport of isolated particles, we explore here a broader range of simulation parameters. In particular, we go beyond the lubrication approximation and investigate systems for which $H/L \simeq 0.02$ up to $H/L \simeq 0.1$. We report the measured long-time diffusion coefficients $D_{eff}$ and mean drift $V_{eff}$ in figure 9, along with representative plots of the density profile in the frame of reference where the channel wall is fixed. We find that as the width $H$ increases, a $y$ dependence of the stationary density emerges (see figure 9e). In fact, diffusion across the channel width can no longer be considered fast with respect to diffusion along length $x$. This corresponds to the progressive breakdown of the lubrication approximation.

Figure 9. Effective (a) diffusion and (b) advection normalised by wall velocity for an ideal Brownian gas as a function of the Péclet number, for different values of the channel height $H$. Error bars correspond to one standard deviation over 10 independent runs. Theory curves correspond to (2.18) and (2.19). Stationary density profiles in the periodic channel for different values of $H$: (c) $H=4\ell _0$, (d) $H=12\ell _0$, and (e) $H=20\ell _0$. (c,d,e) share the same colour scale, where yellow (purple) regions indicate regions of high (low) density, and are all presented for $v_{wall} = 0.5 \ell _0/\tau _0$, corresponding to $Pe \simeq 16$. Numerical parameters are the same as for figure 2, in particular $L = 200\ell _0$.

Surprisingly, the analytic formulas (2.18) and (2.19) are still in remarkable agreement with simulations, up to $H/L \simeq 0.1$. Slight deviations may be observed for $H/L \simeq 0.1$ (corresponding to $H = 20 \times \ell _0$) on the mean drift, where $V_{eff}/v_{wall} > 0$ even at large $Pe$. This is due to accumulation of particles in the upstream bulge, as they collide and leave a wake of particles, instead of having the time to distribute vertically. As a result, some particles are still in the bulge even at large $Pe$, and are therefore carried, which produces a net mean drift. Interestingly, at very small $H/L \simeq 0.02$ (corresponding to $H = 4 \times \ell _0$), we observe slight deviation from the theory this time at small $Pe$. This is due to the fact that for such systems, the penetration depth in the wall, $\delta h \simeq 0.2 \ell _0$, becomes more and more comparable with the channel height $H$. The effective vertical accessible space $H$ is therefore larger, and the value in (2.18) and (2.19) should be modified appropriately (see § A.5).

Appendix C. Perturbation theory to obtain long-time transport coefficients

Here, we perform the perturbation theory to obtain the long-time transport coefficients $D_{eff}$ and ${V_{eff}}$ of the general Fokker–Planck equation (2.12).

In the following, it will be easier to work in Fourier space, and we therefore define, for any function $f(x,t)$, the Fourier transform

(C1)\begin{equation} \tilde{f}(k,\omega) = \int \text{d} x \,\text{d} t \, {\rm e}^{-{\rm i}(kx + \omega t)}\,f(x,t), \end{equation}

where the $\int$ sign encompasses integration over space and time. Conversely, the reverse Fourier transform is given by

(C2)\begin{equation} f(x,t) = \int \frac{\text{d} k\,\text{d} \omega}{(2{\rm \pi})^2} \,{\rm e}^{{\rm i}(kx + \omega t)}\, \tilde{f}(k,\omega). \end{equation}

Performing a Fourier transform on (2.12) yields

(C3)\begin{equation} {\rm i} \omega\,\tilde{p}(k,\omega) =- D_0 k^2\,\tilde{p}(k,\omega) + 1 - {\rm i} \int \frac{\text{d} k' \,\text{d}\omega' }{(2{\rm \pi})^2}\,k\,\tilde{v}(k', \omega')\,\tilde{p}(k-k',\omega - \omega'), \end{equation}

which can be written, using the notation $1/(D_0 k^2 + \textrm {i} \omega )\equiv \tilde {p}_0(k,\omega )$, as

(C4)\begin{equation} \tilde{p}(k,\omega) = \tilde{p}_0(k,\omega) - {\rm i}\,\tilde{p}_0(k,\omega) \int \frac{\text{d} k' \,\text{d}\omega' }{(2{\rm \pi})^2}\,k\,\tilde{v}(k', \omega')\,\tilde{p}(k-k',\omega - \omega'). \end{equation}

We would like to eventually simplify (C4) in a form where diffusion and advection coefficients can be read easily, as suggested by the Fourier transform of (2.11) that yields

(C5)\begin{equation} \tilde{p}_{eff}(k,\omega) = \frac{1}{D_{eff}\,k^2 + {\rm i} (V_{eff}\,k + \omega) }. \end{equation}

The lubrication approximation enables us to simplify further the effective equation on $\tilde {p}(k,\omega )$. Indeed, we consider that the fluctuations at the channel boundary are a perturbation, with small relative amplitude $\varepsilon = h_0/H$. Hence the last term of (C4), containing the advection $v$, can be considered as a perturbation. For example, in the case of the ideal gas,

(C6)\begin{equation} \tilde{v}(k,\omega) \simeq {\rm i} D_0 k\,\frac{\tilde{h}(k,\omega)}{H} = O(\varepsilon), \end{equation}

where we denote $\tilde {h}(k,\omega )$ the Fourier transform of the non-constant part of $h(x,t)$, i.e. $h(x,t)- H$. We therefore seek a solution to (C4) as an expansion in $\varepsilon$, namely, $\tilde {p}(k,\omega ) = \tilde {p}_0(k,\omega ) + \varepsilon \, \tilde {p}_1(k,\omega ) + \varepsilon ^2\,\tilde {p}_2(k,\omega ) +\cdots$.

Additionally, since we seek the behaviour of the solution at long times, it is natural to calculate the noise-averaged solution $\langle \tilde {p}(k,\omega ) \rangle$ where $\langle \cdot \rangle$ denotes the usual average over realisations of the noise. Note that we can also treat the propagating wave case, which is deterministic, in terms of a random field. Consider that we define $h(x,t) = H + h_0 \cos (k_0 x - \omega _0 t + \theta )$, where $\theta$ is a random phase distributed uniformly on $[0,2{\rm \pi} ]$. Clearly, the value of $\theta$ cannot affect the result at late times as it just fixes the initial configuration of the height at time $t = 0$ when the advection diffusion process starts. This choice of $\theta$ is also equivalent, for instance, to choosing an arbitrary initial time $t_0 = \theta /\omega _0$. We thus define $\langle {\cdot } \rangle$ here as a uniform average over $\theta$ on $[0,2{\rm \pi} ]$. Using this convention, we see immediately that $\langle h \rangle = H$. Here, we consider additionally Gaussian fluctuations with mean 0, i.e. the first two moments of the noise completely specify the problem. Thus $\langle \tilde {v}(k,\omega )\rangle = 0$ and

(C7)\begin{equation} \langle \tilde{v}(k,\omega)\,\tilde{v}(k',\omega') \rangle = S(k,\omega)\,(2{\rm \pi})^2 \delta(k + k')\,\delta(\omega + \omega'), \end{equation}

where $S(k,\omega )$ corresponds to the spectrum of the fluctuating advection, as defined in (2.13). Of course, $\langle \tilde {p}_0(k, \omega ) \rangle = \tilde p_0(k, \omega )$.

We now solve iteratively for $\tilde {p}(k,\omega )$. At first order in $\varepsilon$, we obtain

(C8)\begin{equation} \tilde{p}_1(k,\omega) =- {\rm i}\,\tilde{p}_0(k,\omega) \int \frac{\text{d}\omega' \,\text{d} k'}{(2{\rm \pi})^2}\,k\, \tilde{v}(k', \omega')\,\tilde{p}_0(k-k',\omega - \omega') \end{equation}

and $\langle \tilde {p}_1(k,\omega ) \rangle = 0$. At second order in $\varepsilon$, we have

(C9) \begin{align} \tilde{p}_2(k,\omega) &= (- {\rm i})^2\,\tilde{p}_0(k,\omega) \int \frac{\text{d} k' \,\text{d}\omega' }{(2{\rm \pi})^2}\,k \left[ \tilde{v}(k', \omega') \sim p_0(k-k',\omega - \omega') \vphantom{\frac{\text{d} k'' \,\text{d}\omega'' }{(2{\rm \pi})^2}}\right.\nonumber\\ & \quad \left. {}\times \int \frac{\text{d} k'' \,\text{d}\omega'' }{(2{\rm \pi})^2}\,(k-k')\,\tilde{v}(k'', \omega'')\, \sim p_0(k-k'-k'',\omega - \omega'- \omega'')\right], \end{align}

and averaging over noise gives

(C10)\begin{equation} \langle \tilde{p}_2(k,\omega) \rangle = \tilde{p}_0(k,\omega) \left( \int \frac{\text{d} k' \,\text{d}\omega' }{(2{\rm \pi})^2}\,k (k-k')\,\tilde p_0(k-k',\omega - \omega')\,\frac{S(k',\omega')}{h_0^2} \right) \tilde{p}_0(k,\omega). \end{equation}

We observe that we can define $\varSigma (k,\omega )$ such that

(C11)\begin{equation} \varSigma(k,\omega) = \int \frac{ \text{d} k' \,\text{d}\omega'}{(2{\rm \pi})^2}\,\frac{k (k-k') }{D_0 (k-k')^2 + {\rm i} (\omega - \omega')}\,\frac{S(k',\omega')}{h_0^2} \end{equation}

and $\langle \tilde {p}_2(k,\omega ) \rangle = \tilde {p}_0(k,\omega )\,\varSigma (k,\omega )\, \tilde {p}_0(k,\omega )$.

Pursuing the derivation, one can show that the full solution is a geometric series $\langle \tilde {p}\rangle = \tilde {p}_0 + \varepsilon ^2\,\tilde {p}_0 \varSigma \tilde {p}_0 + \varepsilon ^4\, \tilde {p}_0 \varSigma \tilde {p}_0 \varSigma \tilde {p}_0 +\cdots$ that can be re-summed to obtain

(C12)\begin{equation} \langle \tilde{p}(k,\omega) \rangle =\frac{\tilde p_0(k,\omega)}{1-\varepsilon^2\,\varSigma(k,\omega)\,\tilde p_0(k,\omega)}=\frac{1}{D_0 k^2 +\omega^2-\varepsilon^2\,\varSigma(k,\omega)}. \end{equation}

Additional steps may be found in the supplementary information of Marbach et al. (Reference Marbach, Dean and Bocquet2018).

C.1. Long-time results

From the target long-time expression (C5) and from the re-summed propagator (C12), one can read off the effective long-time diffusion coefficient as

(C13)\begin{equation} D_{eff} = D_0 - \frac{\varepsilon^2}{2} \lim_{k\rightarrow0} \partial_{kk} \varSigma(k,0), \end{equation}

and the mean drift as

(C14)\begin{equation} V_{eff} = {\rm i} \varepsilon^2\lim_{k\rightarrow0} \partial_{k} \varSigma(k,0). \end{equation}

Injecting $\varSigma (k,\omega )$ from (C11) into (C13) and (C14) (and dropping the $'$ in the integrals for simplicity), we obtain

(C15)\begin{equation} D_{eff} = D_0 - D_0 \varepsilon^2 \int \frac{\text{d} k \,\text{d}\omega }{D_0(2{\rm \pi})^2}\,\frac{D_0 k^2 + {\rm i}\omega}{(D_0 k^2 - {\rm i}\omega)^2}\,S(k,\omega) \end{equation}


(C16)\begin{equation} V_{eff} =-{\rm i}\varepsilon^2 \int \frac{\text{d} k \,\text{d}\omega }{(2{\rm \pi})^2}\,\frac{k}{D_0 k^2 - {\rm i}\omega}\, S(k,\omega) . \end{equation}

Now, assuming that the problem is translationally invariant in space and time (which is reasonable considering that the channel is assumed to be infinitely long, and that such an assumption still allows one to model many different situations), this implies that the height correlations satisfy $\langle h(x,t)\,h(x',t') \rangle = C(|x-x'|,|t-t'|)$, hence $S(k,\omega ) = S(-k,-\omega )$. In addition, as the correlation function $C$ is real, the conjugate is $S^*(k,\omega ) = S(-k,-\omega )$. Using time and space invariance, we obtain that also $S$ is real. Using the fact that $S$ is even with respect to both its variables, and plugging back the expression of $\varepsilon = h_0/H$, (C15) simplifies to (2.14), and (C16) to (2.15), of the main text.

Appendix D. Additional data for a fluid of interacting particles

We present here additional data for interacting particles. We inspect various values of the interaction strength $\alpha$. We present in figure 10 the results of the normalised effective diffusion $D_{eff}/\overline{D}_0(\rho _0)$ and the normalised mean drift ${V_{eff}} /v_{wall}$. We find very similar results whether probing increasing interaction strength $\alpha$ or probing increasing particle density $\rho _0$ (see figure 3). In fact, increasing $\alpha$ also increases the energy scale $E_0(\alpha,\rho _0)$ contained in a volume $d_c^d$ (with $d$ the dimension of space). Hence, since the theory depends mainly on the value of $E_0= ({\rm \pi} /6) \alpha d^4_c \rho_0$, similar effects are observed naturally when increasing $\alpha$ or $\rho_0$.

Figure 10. Effective (a) diffusion and (b) mean drift for a compressible fluid of soft-core interacting Brownian particles with varying interaction strength $\alpha$. Numerical parameters used here are similar to those in figure 3, and $\rho _0 =1 \ell _0^{-2}$. Error bars correspond to one standard deviation over 10 independent runs. Theory curves for the ideal gas and incompressible fluid are the same as for figures 2(a,b). Theory curves for the fluid of interacting particles correspond to (4.19) and (4.20).


These authors contributed equally to the supervision of the work.


Alexandre, A., Mangeat, M., Guérin, T. & Dean, D.S. 2022 How stickiness can speed up diffusion in confined systems. Phys. Rev. Lett. 128 (21), 210601.CrossRefGoogle Scholar
Alim, K., Amselem, G., Peaudecerf, F., Brenner, M.P. & Pringle, A. 2013 Random network peristalsis in Physarum polycephalum organizes fluid flows across an individual. Proc. Natl Acad. Sci. USA 110 (33), 1330613311.CrossRefGoogle Scholar
Antonov, A.P., Ryabov, A. & Maass, P. 2021 Driven transport of soft Brownian particles through pore-like structures: effective size method. J. Chem. Phys. 155 (18), 184102.CrossRefGoogle Scholar
Arango-Restrepo, A. & Rubi, J.M. 2020 Entropic transport in a crowded medium. J. Chem. Phys. 153 (3), 034108.CrossRefGoogle Scholar
Aris, R. 1956 On the dispersion of a solute in a fluid flowing through a tube. Proc. R. Soc. Lond. A 235 (1200), 6777.Google Scholar
Ben Dor, Y., Ro, S., Kafri, Y., Kardar, M. & Tailleur, J. 2022 Disordered boundaries destroy bulk phase separation in scalar active matter. Phys. Rev. E 105 (4), 044603.CrossRefGoogle Scholar
Bénichou, O., Illien, P., Oshanin, G. & Voituriez, R. 2013 Fluctuations and correlations of a driven tracer in a hard-core lattice gas. Phys. Rev. E 87 (3), 032164.CrossRefGoogle Scholar
Bénichou, O., Illien, P., Oshanin, G., Sarracino, A. & Voituriez, R. 2018 Tracer diffusion in crowded narrow channels. J. Phys.: Condens. Matter 30 (44), 443001.Google Scholar
Bertini, L., De Sole, A., Gabrielli, D., Jona-Lasinio, G. & Landim, C. 2015 Macroscopic fluctuation theory. Rev. Mod. Phys. 87, 593636.CrossRefGoogle Scholar
Bhattacharjee, T. & Datta, S.S. 2019 Bacterial hopping and trapping in porous media. Nat. Commun. 10 (1), 19.CrossRefGoogle ScholarPubMed
Burada, P.S., Schmid, G., Reguera, D., Rubi, J.M. & Hänggi, P. 2007 Biased diffusion in confined media: test of the Fick–Jacobs approximation and validity criteria. Phys. Rev. E 75 (5), 051111.CrossRefGoogle Scholar
Cao, W., Wang, J. & Ma, M. 2019 Water diffusion in wiggling graphene membranes. J. Phys. Chem. Lett. 10 (22), 72517258.CrossRefGoogle Scholar
Chaikin, P.M., Lubensky, T.C. & Witten, T.A. 1995 Principles of Condensed Matter Physics, vol. 10. Cambridge University Press.CrossRefGoogle Scholar
Chakrabarti, B. & Saintillan, D. 2020 Shear-induced dispersion in peristaltic flow. Phys. Fluids 32 (11), 113102.CrossRefGoogle Scholar
Chávez, Y., Chacón-Acosta, G. & Dagdug, L. 2018 Effects of curved midline and varying width on the description of the effective diffusivity of Brownian particles. J. Phys.: Condens. Matter 30 (19), 194001.Google ScholarPubMed
Codutti, A., Cremer, J. & Alim, K. 2022 Changing flows balance nutrient absorption and bacterial growth along the gut. bioRxiv.CrossRefGoogle Scholar
Coquinot, B., Bocquet, L., Kavokine, N. 2023 Quantum feedback at the solid–liquid interface: flow-induced electronic current and its negative contribution to friction. Phys. Rev. X 13 (1), 011019.Google Scholar
Dagdug, L., García-Chung, A.A. & Chacón-Acosta, G. 2016 On the description of Brownian particles in confinement on a non-Cartesian coordinates basis. J. Chem. Phys. 145 (7), 074105.CrossRefGoogle ScholarPubMed
Dean, D.S. 1996 Langevin equation for the density of a system of interacting Langevin processes. J. Phys. A: Math. Gen. 29 (24), L613.CrossRefGoogle Scholar
Dean, D.S. & Lefèvre, A. 2004 Self-diffusion in a system of interacting Langevin particles. Phys. Rev. E 69 (6), 061111.CrossRefGoogle Scholar
Démery, V., Bénichou, O. & Jacquin, H. 2014 Generalized Langevin equations for a driven tracer in dense soft colloids: construction and applications. New J. Phys. 16 (5), 053032.CrossRefGoogle Scholar
Dong, J., Yang, Y. & Zhu, Y. 2022 Recent advances in the understanding of alveolar flow. Biomicrofluidics 16 (2), 021502.CrossRefGoogle ScholarPubMed
Dünweg, B. & Kremer, K. 1991 Microscopic verification of dynamic scaling in dilute polymer solutions: a molecular-dynamics simulation. Phys. Rev. Lett. 66 (23), 2996.CrossRefGoogle ScholarPubMed
Dünweg, B. & Kremer, K. 1993 Molecular dynamics simulation of a polymer chain in solution. J. Chem. Phys. 99 (9), 69836997.CrossRefGoogle Scholar
Evans, J.D., Krause, S. & Feringa, B.L. 2021 Cooperative and synchronized rotation in motorized porous frameworks: impact on local and global transport properties of confined fluids. Faraday Discuss. 225, 286300.CrossRefGoogle ScholarPubMed
Fruchart, M., Scheibner, C. & Vitelli, V. 2023 Odd viscosity and odd elasticity. Annu. Rev. Condens. Matter Phys. 14, 471510.CrossRefGoogle Scholar
Guérin, T. & Dean, D.S. 2015 Kubo formulas for dispersion in heterogeneous periodic nonequilibrium systems. Phys. Rev. E 92 (6), 062103.CrossRefGoogle ScholarPubMed
Haldoupis, E., Watanabe, T., Nair, S. & Sholl, D.S. 2012 Quantifying large effects of framework flexibility on diffusion in MOFs: CH$_4$ and CO$_2$ in ZIF-8. ChemPhysChem 13 (15), 34493452.CrossRefGoogle ScholarPubMed
Hansen, J.-P. & McDonald, I.R. 2013 Theory of Simple Liquids: With Applications to Soft Matter. Academic Press.Google Scholar
Hargus, C., Klymko, K., Epstein, J.M. & Mandadapu, K.K. 2020 Time reversal symmetry breaking and odd viscosity in active fluids: Green–Kubo and NEMD results. J. Chem. Phys. 152 (20), 201102.CrossRefGoogle ScholarPubMed
Hess, W. & Klein, R. 1983 Generalized hydrodynamics of systems of Brownian particles. Adv. Phys. 32 (2), 173283.CrossRefGoogle Scholar
Hydon, P.E. & Pedley, T.J. 1993 Axial dispersion in a channel with oscillating walls. J. Fluid Mech. 249, 535555.CrossRefGoogle Scholar
Illien, P., Bénichou, O., Oshanin, G., Sarracino, A. & Voituriez, R. 2018 Nonequilibrium fluctuations and enhanced diffusion of a driven particle in a dense environment. Phys. Rev. Lett. 120 (20), 200606.CrossRefGoogle Scholar
Jacobs, M.H. 1935 Diffusion Processes. Springer.Google Scholar
Kalinay, P. 2020 Taylor dispersion in Poiseuille flow in three-dimensional tubes of varying diameter. Phys. Rev. E 102 (4), 042606.CrossRefGoogle ScholarPubMed
Kalinay, P. & Percus, J.K. 2006 Corrections to the Fick–Jacobs equation. Phys. Rev. E 74 (4), 041203.CrossRefGoogle Scholar
Kalz, E., Vuijk, H.D., Abdoli, I., Sommer, J.-U., Löwen, H. & Sharma, A. 2022 Collisions enhance self-diffusion in odd-diffusive systems. Phys. Rev. Lett. 129 (9), 090601.CrossRefGoogle ScholarPubMed
Kavokine, N., Bocquet, M.-L. & Bocquet, L. 2022 Fluctuation-induced quantum friction in nanoscale water flows. Nature 602 (7895), 8490.CrossRefGoogle ScholarPubMed
Kavokine, N., Netz, R.R. & Bocquet, L. 2021 Fluids at the nanoscale: from continuum to subcontinuum transport. Annu. Rev. Fluid Mech. 53, 377410.CrossRefGoogle Scholar
Kawasaki, K. 1998 Microscopic analyses of the dynamical density functional equation of dense fluids. J. Stat. Phys. 93, 527546.CrossRefGoogle Scholar
Lahtinen, J.M., Hjelt, T., Ala-Nissila, T. & Chvoj, Z. 2001 Diffusion of hard disks and rodlike molecules on surfaces. Phys. Rev. E 64 (2), 021204.CrossRefGoogle ScholarPubMed
Lappala, A., Zaccone, A. & Terentjev, E.M. 2013 Ratcheted diffusion transport through crowded nanochannels. Sci. Rep. 3 (1), 3103.CrossRefGoogle ScholarPubMed
Lekkerkerker, H.N.W. & Dhont, J.K.G. 1984 On the calculation of the self-diffusion coefficient of interacting Brownian particles. J. Chem. Phys. 80 (11), 57905792.CrossRefGoogle Scholar
Leroy, F., Rousseau, B. & Fuchs, A.H. 2004 Self-diffusion of n-alkanes in silicalite using molecular dynamics simulation: a comparison between rigid and flexible frameworks. Phys. Chem. Chem. Phys. 6 (4), 775783.CrossRefGoogle Scholar
Lowen, H. & Szamel, G. 1993 Long-time self-diffusion coefficient in colloidal suspensions: theory versus simulation. J. Phys.: Condens. Matter 5 (15), 2295.Google Scholar
Ma, M., Grey, F., Shen, L., Urbakh, M., Wu, S., Liu, J.Z., Liu, Y. & Zheng, Q. 2015 Water transport inside carbon nanotubes mediated by phonon-induced oscillating friction. Nat. Nanotechnol. 10 (8), 692695.CrossRefGoogle ScholarPubMed
Ma, M., Tocci, G., Michaelides, A. & Aeppli, G. 2016 Fast diffusion of water nanodroplets on graphene. Nat. Mater. 15 (1), 6671.CrossRefGoogle Scholar
Malgaretti, P., Puertas, A.M. & Pagonabarraga, I. 2022 Active microrheology in corrugated channels: comparison of thermal and colloidal baths. J. Colloid Interface Sci. 608, 26942702.CrossRefGoogle ScholarPubMed
Mangeat, M., Guérin, T. & Dean, D.S. 2017 Dispersion in two dimensional channels – the Fick–Jacobs approximation revisited. J. Stat. Mech.: Theory Exp. 2017 (12), 123205.CrossRefGoogle Scholar
Marbach, S. & Alim, K. 2019 Active control of dispersion within a channel with flow and pulsating walls. Phys. Rev. Fluids 4 (11), 114202.CrossRefGoogle Scholar
Marbach, S. & Bocquet, L. 2019 Osmosis, from molecular insights to large-scale applications. Chem. Soc. Rev. 48 (11), 31023144.CrossRefGoogle Scholar
Marbach, S., Dean, D.S. & Bocquet, L. 2018 Transport and dispersion across wiggling nanopores. Nat. Phys. 14 (11), 11081113.CrossRefGoogle Scholar
Marconi, U.M.B. & Tarazona, P. 1999 Dynamic density functional theory of fluids. J. Chem. Phys. 110 (16), 80328044.CrossRefGoogle Scholar
Masri, R., Puelz, C. & Riviere, B. 2021 A reduced model for solute transport in compliant blood vessels with arbitrary axial velocity profile. Intl J. Heat Mass Transfer 176, 121379.CrossRefGoogle Scholar
Mercer, G.N. & Roberts, A.J. 1990 A centre manifold description of contaminant dispersion in channels with varying flow properties. SIAM J. Appl. Math. 50 (6), 15471565.CrossRefGoogle Scholar
Mercer, G.N. & Roberts, A.J. 1994 A complete model of shear dispersion in pipes. Jpn. J. Ind. Appl. Math. 11, 499521.CrossRefGoogle Scholar
Noh, Y. & Aluru, N.R. 2021 Phonon-fluid coupling enhanced water desalination in flexible two-dimensional porous membranes. Nano Lett. 22 (1), 419425.CrossRefGoogle ScholarPubMed
Obliger, A., Bousige, C., Coasne, B. & Leyssale, J.-M. 2023 Development of atomistic kerogen models and their applications for gas adsorption and diffusion: a mini-review. Energy Fuels 37 (3), 16781698.CrossRefGoogle Scholar
Obliger, A., Valdenaire, P.-L., Ulm, F.-J., Pellenq, R.J.-M. & Leyssale, J.-M. 2019 Methane diffusion in a flexible kerogen matrix. J. Phys. Chem. B 123 (26), 56355640.CrossRefGoogle Scholar
Pàmies, J.C., Cacciuto, A. & Frenkel, D. 2009 Phase diagram of Hertzian spheres. J. Chem. Phys. 131 (4), 044514.CrossRefGoogle ScholarPubMed
Pireddu, G., Pazzona, F.G., Demontis, P. & Załuska-Kotur, M.A. 2019 Scaling-up simulations of diffusion in microporous materials. J. Chem. Theory Comput. 15 (12), 69316943.CrossRefGoogle ScholarPubMed
Puertas, A.M., Malgaretti, P. & Pagonabarraga, I. 2018 Active microrheology in corrugated channels. J. Chem. Phys. 149 (17), 174908.CrossRefGoogle ScholarPubMed
Reguera, D. & Rubí, J.M. 2001 Kinetic equations for diffusion in the presence of entropic barriers. Phys. Rev. E 64 (6), 061106.CrossRefGoogle ScholarPubMed
Reimann, P., Van den Broeck, C., Linke, H., Hänggi, P., Rubi, J.M. & Pérez-Madrid, A. 2001 Giant acceleration of free diffusion by use of tilted periodic potentials. Phys. Rev. Lett. 87 (1), 010602.CrossRefGoogle ScholarPubMed
Riefler, W., Schmid, G., Burada, P.S. & Hänggi, P. 2010 Entropic transport of finite size particles. J. Phys.: Condens. Matter 22 (45), 454109.Google ScholarPubMed
Rizkallah, P., Sarracino, A., Bénichou, O. & Illien, P. 2022 Microscopic theory for the diffusion of an active particle in a crowded environment. Phys. Rev. Lett. 128 (3), 038001.CrossRefGoogle Scholar
Rubi, J.M. 2019 Entropic diffusion in confined soft-matter and biological systems. EPL 127 (1), 10001.