Active Brownian particles near straight or curved walls: Pressure and boundary layers
Ayhan Duzgun Jonathan V. Selinger
cond-mat.soft





Liquid Crystal Institute, Kent State University, Kent, Ohio 44242, USA




Unlike equilibrium systems, active matter is not governed by the conventional laws of thermodynamics. Through a series of analytic calculations and Langevin dynamics simulations, we explore how systems cross over from equilibrium to active behavior as the activity is increased. In particular, we calculate the profiles of density and orientational order near straight or circular walls, and show the characteristic width of the boundary layers. We find a simple relationship between the enhancements of density and pressure near a wall. Based on these results, we determine how the pressure depends on wall curvature, and hence make approximate analytic predictions for the motion of curved tracers, as well as the rectification of active particles around small openings in confined geometries.   


Active Brownian particles near straight or curved walls: Pressure and boundary layers
    Jonathan V. Selinger
    December 27, 2017
=====================================================================================





§ INTRODUCTION


The statistical properties of active particles are quite different from the analogous properties of passive particles <cit.>.  For example, by the conventional laws of thermodynamics, equilibrium Brownian motion cannot perform any useful work.  By contrast, active Brownian motion has been shown to power microscopic gears, thus performing mechanical work <cit.>.  Similarly, active Brownian particles drive the motion of curved tracers <cit.>, and induce flexible membranes to fold <cit.>.  Parallel plates immersed in a bath of active Brownian particles experience an attractive depletion force, analogous to the Casimir effect, unlike plates in a bath of equilibrium Brownian particles <cit.>.  A bath of active particles exerts an active pressure on the walls, which is not a state function of the fluid but rather depends on properties of the walls <cit.>, particularly on the wall curvature <cit.>.  Many of these phenomena arise from a distinctive feature of active systems that self-propelled particles accumulate at walls <cit.>.

The purpose of this article is to explore how systems cross over from equilibrium to active behavior, as a function of particle activity, through several example calculations.  First, in Secs. II and III, we review the derivation of equations of motion for order parameter fields in a system of active particles, and use these equations to calculate the steady-state solutions in several specific geometries:  near a straight wall, between two straight walls, inside and outside a circular wall.  These calculations show the formation of boundary layers near the walls, which are characterized by enhancements in the density and polar order.  We note that similar calculations have been done previously <cit.>; here we review them and verify them through simulation, so that we can apply them to calculate forces and densities in specific geometries in the following sections.  (Most of these calculations are done by truncating a series of equations for moments of the density distribution, but the Appendix shows that some results apply even without this truncation.)

In Sec. IV, we find a simple relationship between the density enhancement in the boundary layer and the active pressure on a wall.  Through this relationship, we determine the pressure on a straight wall, as well as inside and outside a circular wall.  We combine these calculations into a single curvature-dependent active pressure.  This result is consistent with previous results by other investigators, but obtained through a different method.

In Sec. V, we apply these findings to specific geometries that demonstrate important differences between equilibrium and active systems.  For a pair of parallel plates in a bath of active Brownian particles, we use boundary-layer considerations to estimate the Casimir-like depletion force between the plates.  For a curved tracer in an active bath, we find the net force resulting from the different pressures on the inner and outer surfaces.  Finally, for a circular particle corral with just a small opening, we show that particle activity leads to a difference in densities between inside and outside.



§ THEORETICAL FORMALISM


We begin by reviewing the derivation of equations of motion for a system of active Brownian particles, beginning with Langevin dynamics for individual particles and leading up to equations for order parameter fields that can be solved in arbitrary geometries.

We consider a system of active, non-interacting Brownian particles in two dimensions (2D).  We suppose that each particle has a position r=(x,y) and an orientation n̂=(cosθ,sinθ) for its self-propulsive force.  Apart from this self-propulsive force, the particles are isotropic.  To describe the time evolution of the position and orientation, we use Langevin stochastic dynamics in the overdamped limit, which gives



    dr/d t=β D_t [Fn̂-∇U(r)+f(t)],
       dθ/d t=β D_r [g(t)].


Here, D_t and D_r are the translational and rotational diffusion constants, respectively, and β=1/(k_B T) is the inverse temperature.  The coefficient F is the self-propulsive force of each particle, and the product v_0 =β D_t F is the self-propulsive velocity.  The function U(r) is the potential energy, so that -∇U(r) is the force derived from the potential.  The final terms f(t) and g(t) are the random force and torque acting on each particle.  These stochastic contributions satisfy Gaussian white noise statistics, such that

    ⟨ f_i(t)⟩=0,  ⟨ g(t)⟩=0,
       ⟨ f_i(t)f_j(t')⟩=2δ_ijδ(t-t')/β^2 D_t,  ⟨ g(t)g(t')⟩=2δ(t-t')/β^2 D_r,
       ⟨ f_i(t)g(t')⟩=0.


To characterize a large ensemble of active Brownian particles, we use the probability distribution function P(r,θ,t), which gives the probability of finding a particle at position r with orientation θ at time t.  This distribution function evolves in time following the Smoluchowski equation

    ∂ P/∂ t=   -∂/∂ r_i[⟨Δ r_i⟩/Δ tP]
    -∂/∂θ[⟨Δθ⟩/Δ tP]
       +1/2∂^2/∂ r_i ∂ r_j[⟨Δ r_i Δ r_j ⟩/Δ tP]
    +1/2∂^2/∂θ^2[⟨Δθ^2⟩/Δ tP]
       +∂^2/∂ r_i ∂θ[⟨Δ r_i Δθ⟩/Δ tP].

Here, the quantities in angle brackets are averages calculated over all particles at position r with orientation θ during a small time interval from t to t+Δ t.  In our system, direct integration of the Langevin equations gives

    ⟨Δ r_i⟩=(v_0 n_i-β D_t∂_i U)Δ t,  ⟨Δθ⟩=0,
       ⟨Δ r_iΔ r_j⟩=2D_t δ_ijΔ t,  ⟨Δθ^2⟩=2D_r Δ t,
       ⟨Δ r_iΔθ⟩=0.

With those expectation values, the Smoluchowski equation becomes

    ∂ P/∂ t=   -∇·[(v_0n̂-β D_t ∇U)P + D_t ∇P]
    +D_r ∂^2 P/∂θ^2
    
    =   -∇·J+D_r∂^2 P/∂θ^2,

where J(r,θ,t)=(v_0n̂-β D_t ∇U)P + D_t ∇P is the current density of particles at position r with orientation θ at time t.

As a simplification, instead of considering the full distribution P(r,θ,t) as a function of θ, we calculate orientational moments of the distribution


    ρ(r,t)=∫_0^2π P(r,θ,t)dθ,
       M_i(r,t)=∫_0^2π P(r,θ,t)n_i dθ,
       Q_ij(r,t)=∫_0^2π P(r,θ,t)(2n_i n_j -δ_ij)dθ,


with n_1=cosθ and n_2=sinθ.  The zero-th moment ρ(r,t) is the total density of particles, integrated over all orientations, as a function of positition and time.  The higher moments (normalized by ρ) give the orientational order parameters as functions of position and time.  In particular, the vector M(r,t)/ρ(r,t) is the polar order parameter, and the tensor Q(r,t)/ρ(r,t) is the nematic order parameter.  By integrating over the Smoluchowski equation (<ref>), we obtain equations of motion for the moments



    ∂ρ/∂ t=   -∂_i [v_0 M_i -β D_t (∂_i U)ρ -D_t∂_i ρ]
    
    =   -∂_i J^(0)_i,
    ∂ M_j/∂ t=   -∂_i [1/2v_0(ρδ_ij+Q_ij)-β D_t (∂_i U) M_j -D_t ∂_i M_j ]
           -D_r M_j
    
    =   -∂_i J^(1)_ij-D_r M_j.


Here, the moments of current density are defined as J^(0)_i(r,t)=∫_0^2π J_i(r,θ,t)dθ and J^(1)_ij(r,t)=∫_0^2π J_i(r,θ,t)n_j dθ.

In principle, there is an infinite series of equations of motion for the moments, with the equation for the dipole moment M_i depending on the quadrupole moment Q_ij, the equation for the quadrupole moment depending on the octupole moment, and so forth.  As an approximation, we truncate the series by assuming that the quadrupole moment Q_ij=0.  With that approximation, Eqs. (<ref>) provide a closed set of equations for ρ(r,t) and M_i(r,t), which can be solved to find the distribution of active Brownian particles in any geometry.



§ SOLUTION IN SIMPLE GEOMETRIES


In this section, we find steady-state solutions of the Smoluchowski moment equations (<ref>) in several specific geometries with straight or curved walls.  We consider hard walls, so that the potential energy can be written as

    U(r)=
    
    0,    outside wall, 
    ∞,    inside wall.

In the free space outside the wall where U(r)=0, in the steady state, the Smoluchowski moment equations simplify to



    0=-v_0∇·M+D_t ∇^2 ρ,
       0=-1/2v_0 ∇ρ + D_t ∇^2 M - D_r M.


These equations can be combined to give

    ∇^4 ρ = (D_r/D_t+v_0^2/2D_t^2)∇^2 ρ = ξ^-2∇^2 ρ.

Hence, we see that the equation has a natural length scale ξ, which can be written as

    ξ =√(D_t/D_r)/√(1+v_0^2/(2 D_r D_t))=a/√(1+1/2Pe^2).

In the numerator, the ratio of translational and rotational diffusion constants gives the length scale a=√(D_t/D_r), which is typically of the same order as the particle diameter.  The denominator is expressed in terms of the Peclet number Pe=v_0/√(D_t D_r), which is a dimensionless ratio that characterizes the particle activity.

The differential equations must be solved with boundary conditions expressing the constraint that particles cannot enter the hard wall.  These boundary conditions can be written in terms of the moments of current density as



    0=N̂·J^(0)=N̂·(v_0 M-D_t ∇ρ),
       0=N̂·J^(1)=N̂·(1/2v_0 ρI-D_t ∇M),


where N̂ is the local normal to the wall.  This system of equations can be solved exactly in several cases.



 §.§ Particles near an infinite straight wall


Consider an infinite, straight wall along the y-axis, so that the region x<0 is free space, and x>0 is excluded.  By the symmetry of this problem, we assume that ρ and M_x are functions of x only, and M_y=0.  Using these assumptions, we solve the differential equations (<ref>) with the boundary conditions (<ref>), and the additional boundary condition that ρ(x)→ρ _bulk as x→ -∞.  The solution is



    ρ(x)=ρ_bulk(1+v_0^2/2 D_r D_te^x/ξ),
       M(x)=ρ_bulkv_0/2D_r ξe^x/ξ.


This solution is plotted in Fig. <ref> for three sample sets of parameters.  The density ρ(x) is enhanced in a boundary layer of thickness ξ near the wall, and decays exponentially to ρ_bulk.  The maximum density occurs right at the wall, where

    ρ_wall=ρ_bulk(1+v_0^2/2 D_r D_t).

In the Appendix, we will show that this result is exact, and does not depend on the truncation of moments.

The first moment M_x(x) is nonzero in the same boundary layer, and decays exponentially to zero.  Because M_x(x) is positive, we can see that particles accumulate at the wall with their orientations pointing into the wall, but are unable to enter the wall.  The orientational order parameter M_x(x)/ρ(x) is greatest at the wall, and decays exponentially into the bulk, which is isotropic.



As a numerical check on the calculation, we perform simulations of the Langevin equations (<ref>) with a boundary at x=0, and histograms of the density ρ(x) and first moment M_x(x) are plotted in Fig. <ref>.  The numerical results agree very well with the analytic predictions.  



 §.§ Particles between two walls




Now consider a system with two infinite, parallel walls at x=± L.  We solve the differential equations (<ref>) with the boundary conditions (<ref>) on both walls.  The solution is



    ρ(x)=ρ_1[1+v_0^2/2 D_r D_tcosh(x/ξ)/cosh(L/ξ)],
       M_x(x)=ρ_1 v_0/2 D_rξsinh(x/ξ)/cosh(L/ξ),


with an overall coefficient ρ_1.  This solution is plotted in Fig. <ref>.  Here, the density ρ(x) has boundary layers of thickness ξ near both walls, and the first moment M_x(x) points into each of the walls.

Instead of a constraint on the density far from the wall, we now have a constraint on the integrated number of particles in the system, which can be written as ∫_-L^L ρ(x)dx=2Lρ̅, where ρ̅ is the average density.  Hence, the overall coefficient is

    ρ_1 = ρ̅[1+v_0^2 ξ/2 D_r D_t Ltanh(L/ξ)]^-1.

For large wall separation, with L≫ξ[1+v_0^2 /(2 D_r D_t)], the density in the center is approximately independent of the walls, and we can just write the overall coefficient as ρ_1=ρ_bulk.  However, for smaller wall separation, the density in the center is depleted because the density on the walls is enhanced, as shown in the figure.

We perform numerical simulations of the Langevin equations with boundaries on both side of the domain.  In the simulation algorithm, when a particle crosses a boundary, we place it back in its previous location and update the orientation through the usual rotational diffusion.  The numerical results, plotted in Fig. <ref>, agree very well with the analytic predictions.



 §.§ Particles inside circle


Suppose that active particles are confined inside a hard, circular wall of radius R.  In this case, it is most convenient to work in terms of polar coordinates (r,θ).  By rotational symmetry, we expect that ρ and M_r are functions of r only, and M_θ =0.  We then express the differential equations (<ref>) in terms of polar coordinates.  A general solution for ρ(r) is a linear combination of the modified Bessel functions I_0(r/ξ) and K_0(r/ξ), and the corresponding solution for M_r(r) is a linear combination of I_1(r/ξ) and K_1(r/ξ).  Because the density must be finite at r=0, ρ(r) can have only the I_0(r/ξ) function, and hence M_r(r) can have only I_1(r/ξ).  By putting these functions into the boundary conditions (<ref>), expressed in polar coordinates, we obtain



    ρ(r)=ρ_2[1+v_0^2/D_r D_tI_0 (r/ξ)/A_2],
       M_r(r)=ρ_2 v_0/D_r ξI_1 (r/ξ)/A_2,

with the denominator

    A_2=(1+v_0^2/2D_r D_t)I_2(R/ξ)+(1-v_0^2/2D_r D_t)I_0(R/ξ),


and with an overall coefficient ρ_2.  This solution for density ρ(r) shows a boundary layer of thickness ξ inside the circular wall.  The first moment M_x(x) points outward from the center, into the circular wall.

As in the previous case, the overall coefficient ρ_2 is fixed by the constraint on the integrated number of particles in the system, ∫_0^R2π r ρ(r) dr=π R^2 ρ̅, which gives

    ρ_2 = ρ̅[1+2 v_0^2 ξ I_1 (R/ξ)/D_r D_t R A_2]^-1.

For a large circle, with R≫ξ[1+v_0^2 /(2 D_r D_t)], the density in the center is approximately independent of the wall, and we can write the overall coefficient as ρ_2=ρ_bulk.  However, for smaller radius, the density in the center is depleted because the density at the wall is enhanced.



 §.§ Particles outside circle


As a modification of the previous case, consider active particles that are confined outside a hard, circular wall of radius R.  Again, we work in polar coordinates and express the solution in terms of modified Bessel functions.  Now the density must be finite, with ρ(r)→ρ_bulk, as r→∞.  Hence, ρ(r) can include only the K_0(r/ξ) function, and M_r(r) can include only K_1(r/ξ).  From the boundary conditions on the circular wall, the solution becomes



    ρ(r)=ρ_bulk[1+v_0^2/D_r D_tK_0(r/ξ)/A_3],
       M_r(r)=-ρ_bulk v_0/D_r ξK_1(r/ξ)/A_3,

with the denominator

    A_3=(1+v_0^2/2D_r D_t)K_2(R/ξ)+(1-v_0^2/2D_r D_t)K_0(R/ξ).


In this solution, the density ρ(r) has a boundary layer of thickness ξ outside the circular wall.  The first moment M_x(x) points in toward the center, into the circular wall.



§ PRESSURE ON STRAIGHT OR CURVED WALLS


In the previous section, we calculated the enhancement of density in boundary layers along hard walls of different shapes:  straight, inside a circle, and outside a circle.  In addition to the density enhancement, we would also like to calculate the pressure of active particles against each of these walls.  To calculate the pressure, we use a method based on the theory of Solon et al. <cit.>, and consider a hard wall to be the limiting case of a soft wall.

As a first step, consider an infinite straight wall along the y-axis.  Instead of the hard wall potential energy of Eq. (<ref>), we use the soft potential

    U(x)=
    
    0,    for  x<0, 
    
    U' x,    for  x>0,

where U' is a finite positive constant.  The limit of U'→∞ will then represent a hard wall.  This assumption is similar to Ref. <cit.>, except that we use a linear potential while they used a quadratic potential.

We now solve the steady-state differential equations (<ref>) separately in the regions x<0 and x>0.  In each region, we look for solutions where the density ρ(x) and the first moment M_x(x) vary as e^α x.  The differential equations then give a characteristic equation for α.  In the region x<0, the characteristic equation is

    α^3-(D_r/D_t+v_0^2/2D_t^2)α=0,

and the solutions are α=0 or ±1/ξ.  In the region x>0, the characteristic equation is 

    α^3+2β U' α^2-(D_r/D_t+v_0^2/2D_t^2-β^2 U'^2)α-D_r β U'/D_t=0.

For large U', the solutions are

    α=D_r/β U' D_t+O(U'^-2) or -β U' ±v_0/√(2)D_t+O(U'^-1).

Because ρ cannot diverge as x→±∞, we must eliminate the negative value of α for x<0, and the positive value of α for x>0.  We then have four exponential modes with coefficients to be determined from the boundary conditions.  At the boundary x=0, we require that the density ρ(x), the first moment M_x(x), and the current moments J^(0)_x(x) and J^(1)_xx(x) must all be continuous (keeping in mind that the definition of current moments in Eqs. (<ref>) includes U' terms in the region x>0).  We also require that ρ(x)→ρ_bulk as x→-∞, far from the wall.  Applying these boundary conditions, and assuming that β U' D_t /v_0 ≫ 1, we obtain

    ρ(x)=
    ρ_bulk(1+v_0^2/2D_r D_te^x/ξ),    for  x<0, 
    ρ_bulk(1+v_0^2/2D_r D_t)e^-β U' x,    for  x>0.

The result for x<0 is exactly the same as previously calculated for a hard wall in Eq. (<ref>a).  The result for x>0 shows how the density decreases inside the wall, dominated by the Boltzmann distribution for large U'.

From these results, we can calculate the pressure of the particles on the wall.  As noted in Ref. <cit.>, the force of the particles on the wall is equal and opposite to the force of the wall on the particles.  Hence, the pressure can be calculated as

    p   =∫_0^∞ρ(x)∂ U(x)/∂ xdx=U'∫_0^∞ρ(x)dx
       =ρ_bulk/β(1+v_0^2/2 D_r D_t)
    =ρ_bulk k_B T (1+v_0^2/2 D_r D_t).

This result is consistent with Ref. <cit.>.  In Eq. (<ref>), the first term is the pressure of an ideal gas without activity, p_ideal=ρ_bulk k_B T.  The second term is an enhancement due to the active velocity v_0.  Hence, the active pressure is enhanced over the ideal gas pressure by a factor of (1+v_0^2/(2D_r D_t)).  By comparison, Eq. (<ref>) shows that the density at the wall is enhanced over the bulk density by the same factor ρ_wall=ρ(0)=ρ_bulk(1+v_0^2/(2D_r D_t)).  Hence, the active pressure is simply related to the enhanced density at the wall by p=ρ_wall k_B T.

This relationship between pressure and enhanced density at the wall is quite general.  If we only assume that the wall potential U(x) is large, diverging as x→∞, so that the density inside the wall is dominated by the Boltzmann distribution ρ(x)=ρ_wall e^-β U(x), then the pressure becomes

    p   =∫_0^∞ρ(x)∂ U(x)/∂ xdx=ρ_wall∫_0^∞ e^-β U(x)∂ U(x)/∂ xdx
       =ρ_wall∫_0^U(∞) e^-β UdU = ρ_wall/β=ρ_wall k_B T.


Now we apply the same considerations to active particles inside or outside a circular wall.  For particles inside a circle, the density profile ρ(r) is given by Eq. (<ref>a), with the denominator A_2 defined in Eq. (<ref>c).  Suppose the radius R is large enough that the coefficient ρ_2 can be approximated by ρ_bulk.  The density ρ_wall is just the value of ρ(r) at r=R, and hence the pressure on the wall is

    p=ρ_bulkk_B T[1+v_0^2/D_r D_t/1-v_0^2/2D_r D_t
    +(1+v_0^2/2D_r D_t)I_2 (R/ξ)/I_0 (R/ξ)].

For R≫ξ, we can use the asymptotic expansion I_2(x)/I_0(x)≈ 1-2x^-1+x^-2 to obtain

    p=ρ_bulkk_B T[1+v_0^2/2 D_r D_t+v_0^2/2 D_r^2 ξ R+v_0^4 + D_r D_t v_0^2/4 D_r^3 D_t R^2].

Hence, the pressure of active particles inside a circular wall is increased, compared with the active pressure on a straight wall, by an amount proportional to 1/R.

For particles outside a circular wall, the density profile ρ(r) is given by Eq. (<ref>a), with the denominator A_3 defined in Eq. (<ref>c).  The density ρ_wall is just the value of ρ(r) at r=R, and hence the pressure on the wall is

    p=ρ_bulkk_B T[1+v_0^2/D_r D_t/1-v_0^2/2D_r D_t
    +(1+v_0^2/2D_r D_t)K_2 (R/ξ)/K_0 (R/ξ)].

For R≫ξ, we can use the asymptotic expansion K_2(x)/K_0(x)≈ 1+2x^-1+x^-2 to obtain

    p=ρ_bulkk_B T[1+v_0^2/2 D_r D_t-v_0^2/2 D_r^2 ξ R+v_0^4 + D_r D_t v_0^2/4 D_r^3 D_t R^2].

Hence, the pressure of active particles outside a circular wall is reduced, compared with the active pressure on a straight wall, by an amount proportional to 1/R.

The three cases of straight wall, inside circle, and outside circle can all be combined into the single concept of a curvature-dependent active pressure.  Let us define the curvature κ as κ=0 for a straight wall, κ=+1/R inside a circle, and κ=-1/R outside a circle.  The three equations (<ref>), (<ref>), and (<ref>) can then be combined into the equation

    p(κ)=ρ_bulkk_B T[   (1+v_0^2/2 D_r D_t)+(v_0^2/2 D_r^2 ξ)κ
       +(v_0^4 + D_r D_t v_0^2/4 D_r^3 D_t)κ^2].

This expression can be written more compactly in terms of the Peclet number Pe=v_0/√(D_t D_r) and the particle length scale a=√(D_t/D_r) as

    p(κ)=ρ_bulkk_B T[   (1+1/2Pe^2)
    +1/2Pe^2 √(1+1/2Pe^2) aκ
       +1/4Pe^2(1+Pe^2)a^2 κ^2].

From this result, we can see that the pressure of an active fluid on a wall depends on the shape of the wall through the curvature κ.  Of course, this is not the case for an equilibrium fluid; the pressure of an equilibrium fluid depends only on fluid properties, not on the shape of the wall.

This result for the curvature-dependent active pressure is similar to a previous theoretical result of Fily et al. <cit.>.  They investigated the active pressure inside a convex box with variable positive curvature, and found that the local pressure is proportional to the local curvature of the boundary.  Here we see that the result also applies in a region of negative curvature, such as the outside of a circular wall.



§ APPLICATIONS


In this section we discuss three examples, where we can make predictions for the behavior of active systems based on the concepts of pressure and boundary layers.



 §.§ Depletion force between two plates


Consider two parallel plates inside a bath of active particles, separated by a distance of 2L, as shown in Fig. <ref>.  This problem has been investigated through simulations in Ref. <cit.>.  Outside the plates, far from the ends, this problem is equivalent to the infinite straight wall discussed in Sec. III(A), and the density profile is given by Eq. (<ref>a).  Between the plates, far from the ends, this problem is equivalent to the infinite parallel walls discussed in Sec. III(B), and the density profile is given by Eq. (<ref>a).  Hence, there is a boundary layer with enhanced density on each side of each plate, and thus there is active pressure on each side of each plate.  The question is:  Do all the boundary layers have the same density enhancement?  If the inner boundary layers have a different density enhancement than the outer boundary layers, then the active pressure will either push the plates together or push them apart.



To answer that question, we note that the density profiles (<ref>a) and (<ref>a) have two different overall coefficients.  For the density profile outside the plates, the coefficient is ρ_bulk, which is the density far from the plates.  For the density profile between the plates, the coefficient is ρ_1.  To determine ρ_1, we assume that the region midway between the plates is in contact with the bulk region through the openings at the ends, so that the density midway between the plates is equal to ρ_bulk.  This assumption implies that the density profile between the walls is

    ρ_in(x)=ρ_bulk/1+v_0^2/2 D_r D_tsech(L/ξ)[1+v_0^2/2 D_r D_tcosh(x/ξ)/cosh(L/ξ)].

Hence, the active pressure on the inside surface of each wall is

    p_in=ρ_in^wallk_B T
    =ρ_bulk k_B T/1+v_0^2/2 D_r D_tsech(L/ξ)(1+v_0^2/2 D_r D_t).

By comparison, the active pressure on the outside surface of each wall is

    p_out=ρ_out^wallk_B T
    =ρ_bulk k_B T (1+v_0^2/2 D_r D_t).

Thus, the inside pressure is less than the outside pressure by

    Δ p(L)=p_out-p_in
    =ρ_bulk k_B T (v_0^2/2 D_r D_t)(1+v_0^2/2 D_r D_t)/cosh(L/ξ)+v_0^2/2 D_r D_t.
 
This pressure difference pushes the plates together.  It decays exponentially with L, with the characteristic length scale ξ.  Physically, this force on the plates can be regarded as a depletion force, associated with the reduced boundary layer between the plates compared with outside the plates.  It appears analogous to the Casimir force between conducting plates, but arises from a different mechanism.

We perform Langevin dynamics simulations of a bath of active Brownian particles around two parallel plates, illustrated in Fig. <ref>.  These simulations show that boundary layers form on both sides of both plates, and the outer boundary layers have higher density than the inner boundary layers, as indicated by the red line in the inset.  The relative density of these two boundary layers depends on the separation 2L between the plates.  The main figure presents the density difference, which is proportional to the pressure difference, in comparison with the prediction of Eq. (<ref>).  We can see that the trends are consistent for large separation.  For smaller separation, the prediction overestimates the density difference, perhaps because it is more difficult for the density midway between the plates to become equal with ρ_bulk when the openings between the plates are so small.



 §.§ Force on a curved tracer particle


Consider a curved tracer surrounded by a bath of active particles, as shown in Fig. <ref>.  This type of geometry has been studied through simulations in Refs. <cit.>.  A boundary layer forms on both sides of the tracer, and it experiences active pressure on both sides.  Based on the argument in Sec. IV, the pressure in the inner side is greater than the pressure on the outer side.  As a result, the bath of active particles exerts a net force on the tracer, causing it to move.



To estimate the net force, we use Eq. (<ref>) for the pressure as a function of curvature.  On the inner side, we have the curvature κ=1/(R-s/2), where R is the radius of the midline and s is the thickness of the tracer.  On the outer side, we have κ=-1/(R-s/2).  We assume that ρ_bulk is the same on both sides, because the bulk regions can easily exchange particles.  Hence, for large R and small s, the net pressure becomes

    Δ p=ρ_bulkk_B T v_0^2/D_r^2 ξ R(1+s^2/R^2).


Through Langevin dynamics simulations, we visualize the distribution of active particles around the tracer and calculate the density along the symmetry axis, as indicated by the red line in the Fig. <ref> inset.  This simulation result shows the higher density on the inner side than on the outer side, and hence a higher pressure.  In the main figure, we show the density difference in the simulation in comparison with the prediction from Eq. (<ref>).  This results show a consistent trend, although the prediction is higher by about a factor of 2.  Hence, the approximate argument about a curvature-dependent pressure provides a simple way to understand the net pressure on the tracer.

Experiments have demonstrated that the active motion of swimming bacteria causes an asymmetric gear to rotate <cit.>.  The structure of the asymmetric gear is equivalent to several curved tracer particles linked together, and hence we expect that the argument in this section would also apply to the gear.



 §.§ Corral


As a variation on the curved tracer, consider the active particle corral shown in Fig. <ref>.  Here, the wall is almost a full circle, with only a small opening connecting inside and outside.  The density at the inner boundary is related the the density at the center of the  circle by Eq. (<ref>a), and the density at the outer boundary is related to the bulk density by Eq. (<ref>a).  The question is then:  How are the inside and outside densities related to each other?  For this geometry with only a small opening, it seems reasonable that the density of the inside boundary layer should match the density of the outside boundary layer, ρ_in(R)=ρ_out(R).  This equation gives us a relationship between the density at the center of the circle and the bulk density.  For R≫ξ, that relationship becomes

    ρ_center/ρ_bulk=1-v_0^2/D_r D_tξ/R.

This equation shows that the corral has a reduced density at the center, compared with the bulk density.  Intuitively, this behavior occurs because it is easier for active particles to escape from the corral than to enter it, because of the shape of the opening.  It would not occur for an equilibrium fluid, which would have equal density inside and outside the corral.  Thus, the corral geometry provides one potential mechanism for tunable rectification of active particles, analogous to other mechanisms suggested by Ref. <cit.>.



We carry out Langevin dynamics simulations of a bath of active particles around the corral in Fig. <ref>.  These simulations show that the average density around the center is reduced compared with the bulk exterior, as shown by the red line.  The main figure presents the simulation results for the density ratio as a function of activity, in comparison with the prediction from Eq. (<ref>), and these results are generally consistent.  Hence, the approximate argument about boundary layers provides a way to understand the relative densities in this geometry.



We would like to thank A. Baskaran for helpful discussions.  This work was supported by National Science Foundation Grant No. DMR-1409658. 

*


§ 


The purpose of this appendix is to show that the relationship between the wall density ρ_wall and the bulk density ρ_bulk is exact in one dimension, and does not depend on the truncation of moments.

We begin with the moment equations (<ref>) in steady state,



    0=   -∂ J^(0)_x/∂ x,
    
    0=   -∂ J^(1)_xx/∂ x-D_r M_x .


Integrating Eq. (<ref>a) gives J^(0)_x = constant.  Because no particles are entering or leaving the system at x→±∞, we must have J^(0)_x = 0.  From the definition of J^(0)_x in Eq. (<ref>a), this equation implies

    v_0 M_x - D_t ∂ρ/∂ x = 0.

Integrating once again, we obtain

    v_0 ∫_0^∞ M_x (x) dx = D_t [ρ(∞)-ρ(0)].

By comparison, integrating Eq. (<ref>b) implies

    D_r ∫_0^∞ M_x (x) dx= J^(1)_xx(0)-J^(1)_xx(∞).

At the wall x=0, there is no current, so that J^(1)_xx(0)=0.  Away from the wall, the fluid becomes isotropic, and hence all the moments of the distribution function vanish, except for ρ(∞)=ρ_bulk.  From the definition of J^(1)_xx in Eq. (<ref>b), we obtain J^(1)_xx(∞)=1/2v_0ρ_bulk.  Combining Eqs. (<ref>) and (<ref>) then yields

    D_t/v_0[ρ_bulk-ρ(0)]=-v_0/2D_rρ_bulk,
       ρ_wall=ρ(0)=ρ_bulk(1+v_0^2/2 D_r D_t).

Hence, the active pressure directly follows as

    p=ρ_wallk_B T=ρ_bulkk_B T(1+v_0^2/2 D_r D_t).




