Focusing of Active Particles in a Converging Flow

We consider active particles swimming in a convergent fluid flow in a trapezoid nozzle with no-slip walls. We use mathematical modeling to analyze trajectories of these particles inside the nozzle. By extensive Monte Carlo simulations, we show that trajectories are strongly affected by the background fluid flow and geometry of the nozzle leading to wall accumulation and upstream motion (rheotaxis). In particular, we describe the non-trivial focusing of active rods depending on physical and geometrical parameters. It is also established that the convergent component of the background flow leads to stability of both downstream and upstream swimming at the centerline. The stability of downstream swimming enhances focusing, and the stability of upstream swimming enables rheotaxis in the bulk.


Introduction
Active matter consists of a large number of self-driven agents converting chemical energy, usually stored in the surrounding environment, into mechanical motion [1,2,3]. In the last decade various realizations of active matter have been studied including living self-propelled particles as well as synthetically manufactured ones. Living agents are for example bacteria [4,5], microtubules in biological cells [6,7], spermatozoa [8,9,10] and animals [11,12,13].
These phenomena observed in living active matter can also be achieved using synthetic swimmers, such as self-thermophoretic [55] and self-diffusiophoretic [56,57,58,59] micron sized particles as well as particles set into active motion due to the influence of an external field [60,61,62].
Using simple models we describe the extrusion of a dilute active suspension through a trapezoid nozzle. We analyze the qualitative behavior of trajectories of an individual active particle in the nozzle and study the statistical properties of the particles in the nozzle. The accumulation at walls and rheotaxis are important for understanding how an active suspension is extruded through a nozzle. Wall accumulation may eliminate all possible benefits caused by the activity of the particles in the bulk. Due to rheotaxis active particles may never reach the outlet and leave the nozzle through the inlet, so that properties of the suspension coming out through the outlet will not differ from those of the background fluid.
The specific geometry of the nozzle is also important for our study. The nozzle is a finite domain with two open ends (the inlet and the outlet) and the walls of the nozzle are not parallel but convergent, that is, the distance between walls decreases from the inlet to the outlet. The statistical properties of active suspension (e.g., concentration of active particles) extruded in the infinite channel with parallel straight or periodic walls are well-established, see e.g., [63] and [64], respectively. The finite nozzle size leads to a "proximity effect", i.e., the equilibrium distribution of active particles changes significantly in proximity of both the inlet and the outlet. The fact that the walls are convergent, results in a "focusing effect", i.e., the background flow compared to the pressure driven flow in the straight channel (the Poiseuille flow) has an additional convergent component that turns a particle toward the centerline. Specifically, in this work it is shown that due to this convergent component of the background flow both up-and downstream swimming at the centerline are stable. Stability of the upstream swimming at the centerline is somewhat surprising since from observations in the Poisueille flow it is expected that an active particle turns against the flow only while swimming towards the walls, where the shear rate is higher. This means that we find rheotaxis in the bulk of an active suspension.

Model
To study the dynamics of active particles in a converging flow, two modeling approaches are exploited. In both, an active particle is represented by a rigid rod of length swimming in the xy-plane. In the first -simpler -approach, the rod is a one-dimensional segment which cannot penetrate a wall, whereas in the second -more sophisticated -approach we use the Yukawa segment model [65] to take into account both finite length and width of the rod, as well as a more accurate description of particle-wall steric interaction.
The active particle's center location and its unit orientation vector are denoted by r = (x, y) and p = (cos ϕ, sin ϕ), respectively. The active particles are self-propelled with a velocity directed along their orientation v 0 p. The active particles are confined by a nozzle, see Fig. 1, which is an isosceles trapezoid Ω, placed in the xy-plane so that inlet x = x in and outlet x = x out are bases and the y-axis is the line of symmetry: The nozzle length, the distance between the inlet and the outlet, is denoted by L, i.e., L = |x out − x in |. The width of the outlet and the inlet are denoted by w out and w in , respectively, and their ratio is denoted by k = w out /w in . Furthermore, the active particles are exposed to an external background flow. We approximate the resulting converging background flow due to the trapezoid geometry of the nozzle by where u 0 is a constant coefficient related to the flow rate and α is the slope of walls of the nozzle. Equation (2) is an extension of the Poiseuille flow to channels with convergent walls 1 . Active particles swim in the low Reynolds-number regime. The corresponding overdamped equations of motion for the locations r and orientations p are given by: Here (4) is the Jeffery's equation [20,66,67] for rods with an additional term due to random re-orientation with rotational diffusion coefficient D r ; ζ is an uncorrelated noise with the intensity ζ(t), ζ(t ) = δ(t − t ), e ϕ = (− sin ϕ, cos ϕ). Equation (4) can also be rewritten for the orientation angle ϕ: Here ω = 1 2 , and γ = 1 2 ∂u y ∂x + ∂u x ∂y are local vorticity, vertical expansion (or, equivalently, horizontal compression; similar to Poisson's effect in elasticity) and shear.
The strength of the background flow is quantified by the inverse Stokes number, which is the ratio between the background flow at the center of the inlet and the self-propulsion velocity v 0 . Specifically, where (x in , 0) denotes the location at the center of the inlet.
In the first modeling approach we include the particle wall interaction in the following way: an active particle is not allowed to penetrate the walls of the nozzle. To enforce this, we require that both the front and the back of the particle, r(t) ± ( /2)p, are located inside the nozzle. In numerical simulations of the system (3)-(5) this requirement translates into the following rule: if during numerical integration of (3)-(5) a particle penetrates one of the two walls, then this particle is instantaneously shifted back along the inward normal at the minimal distance, so its front and back are again located inside the nozzle while its orientation is kept fixed.
Unless mentioned otherwise, in this modeling approach we consider a nozzle whose inlet width w in = 0.2 mm and outlet width w out = 0.1 mm are fixed. The following nozzle lengths are considered: L = 0.2 mm, L = 0.5 mm and L = 1.0 mm. The length of the active particles is = 20 µm, they swim with a self-propulsion velocity v 0 = 10 µm s −1 and their rotational diffusion coefficient is given by D r = 0.1 s −1 .
All active particles are initially placed at the inlet, x(0) = x in , with random y-component y(0) and orientation angle ϕ(0). The probability distribution function for initial conditions y(0) and ϕ(0) is given by Ψ ∝ 1 (uniform). The trajectory of an active particle is studied until it leaves the nozzle either through the inlet or the outlet. To gather statistics we use 96,000 trajectories.
We use the second approach to describe the particle-wall interactions and the torque induced by the flow more accurately. For this purpose each rod, representing an active particle, of length , width λ and the corresponding aspect ratio a = /λ is discretized into n r spherical segments with n r = 9a/8 ( x denotes the nearest integer function). The resulting segment distance is also used to discretize the walls of the nozzle into n w segments in the same way. Between the segments of different objects a repulsive Yukawa potential is imposed. The resulting total pair potential is given by where λ is the screening length defining the particle diameter, U 0 is the prefactor of the Yukawa potential and r ij = |r i − r j | is the distance between segment i of a rod and j of the wall of the nozzle, see Fig. 2. Figure 2: Sketch of a discretized active rod (red) of length and width λ which is propelled with a velocity v 0 along its orientation p and is exposed to a converging background flow u BG in the presence of a trapezoid nozzle confinement of length L and with an inlet of size w in and outlet of size w out (blue). To study a system with a packing fraction ρ = 0.1 a channel is attached to the inlet with a non-converging background flow.
The equations of motion (3) and (4) are complemented by the respective derivative of the total potential energy of a rod along with the one-body translational and rotational friction tensors for the rods f T and f R which can be decomposed into parallel f , perpendicular f ⊥ and rotational f R contributions which depend solely on the aspect ratio a [68]. For this approach we measure distances in units of λ, velocities in units of v 0 = F 0 /f (here F 0 is an effective self-propulsion force), and time in units of τ = λf /F 0 . While the width of the outlet w out is varied, the width of the inlet w in as well as the length of the nozzle L is fixed to 100λ in our second approach. Initial conditions are the same as in the first approach. To avoid that a rod and a wall initially intersect each other, the rod is allowed to reorient itself during an equilibration time t e = 10τ while its center of mass is fixed.
Furthermore, we use the second approach to study the impact of a finite density of swimmers. For this approach we initialize N active rods in a channel confinement which is connected to the inlet of the nozzle, see Fig. 2. Inside the channel we assume a regular (non-converging) Poiseuille flow [69]. We restrict our study to a dilute active suspension with a two dimensional packing fraction ρ = 0.1. To maintain this fraction, particles which leave the simulation domain are randomly placed at the inlet of the channel confinement.

Focusing of outlet distribution
Here we characterize the properties of the particles leaving the nozzle at either the outlet or the inlet. Specifically, our objective is to determine whether particles accumulate at the center or at walls when they pass through the outlet or the inlet.
We start with the first modeling approach. Figure 3 shows the spatial distribution of active particles leaving the nozzle at the outlet for various inverse Stokes number σ and three different lengths L of the nozzle, while the width of the inlet and the outlet are fixed. For small  inverse Stokes number σ, the background flow is negligible compared to the self-propulsion velocity. Active particles swim close to the walls and peaks at walls are still clearly visible for σ = 0.5 for all nozzle lengths L, see Fig. 3(a). For σ = 1, the self-propulsion velocity and the background flow are comparable; in this case the histogram shows a single peak at the center of the outlet, see Fig. 3(b). Further increasing the inverse Stokes number from σ = 1 to σ = 9 leads to a broadening of the central peak and then to the formation of two peaks with a well in the center of the outlet, see Fig. 3(c)-(e). Finally, for an even larger inverse Stokes number σ, the self-propulsion velocity is negligible and the histogram becomes close to the one in the passive (no self-propulsion, v 0 = 0) case, see Fig. 3(f). Here the histogram for a nozzle length L = 0.2 mm is uniform except at the edges where it has local peaks due to accumulation at the walls caused by steric interactions.
Histograms for both the y-component and the orientation angle ϕ of the active particles reaching the outlet are depicted in Fig. 4(a)-(c). While active particles leave the nozzle with orientations away from the centerline for small inverse Stokes number, σ = 0.5, they are mostly oriented towards the centerline for larger values of the inverse Stokes number. In Fig. 4(c), one can observe that the histogram is concentrated largely for downstream orientations ϕ ≈ 0 and slightly for upstream orientations ϕ ≈ ±π. These local peaks for ϕ ≈ ±π away from walls are evidence of rheotaxis in the bulk. These peaks are visible for large inverse Stokes numbers only and the corresponding active particles are flushed out of the nozzle with upstream orientations. Due to rotational diffusion and rheotaxis it is possible that an active particle can leave the nozzle through the inlet. We compute the probability of active particles to reach the outlet. This probability, as a function of the inverse Stokes number σ for the three considered nozzle lengths L, is shown in Fig. 5(a), together with selected trajectories, see insets in Fig. 5(a). The figure shows that the probability that an active particle eventually reaches the outlet monotonically grows with the inverse Stokes number σ. Note that a passive particle always leaves the nozzle through the outlet. By comparing the probabilities for different nozzle lengths L it becomes obvious that an active particle is less likely to leave the nozzle through the outlet for longer nozzles. Due to the larger distance L between the inlet and the outlet an active particle spends more time within the nozzle, which makes it more likely to swim upstream by either rotational diffusion or rheotaxis. In Fig. 5(b)-(d) histograms for active particles leaving the nozzle through the inlet are shown. In the case of small inverse Stokes number, σ = 0.5, the majority of active particles leaves the nozzle at the inlet. Specifically, most of them swim upstream due to rheotaxis close to the walls, but some active particles leave the nozzle at the inlet close to the center. These active particles are oriented upstream due to random reorientation. By increasing the inverse Stokes number σ ≥ 1, active particles are no longer able to leave the nozzle at the inlet close to the center. Figure 6: Examples of two trajectories for L = 1 mm and σ = 1.0. The red trajectory starts and ends at the inlet (the endpoint is near the lower wall). The blue trajectory has a zigzag shape with loops close to the walls; the particle that corresponds to the blue trajectory manages to reach the outlet.
Let us now consider specific examples of active particles' trajectories, see Fig. 6. The first trajectory (red) starts and ends at the inlet. Initially the active particle swims downstream and collides with the upper wall due to the torque induced by the background flow. Close to the wall it exhibits rheotactic behavior, but before it reaches the inlet it is expelled towards the center of the nozzle due to rotational diffusion, similar to bacteria that may escape from surfaces due to tumbling [70]. Eventually, the active particle leaves the nozzle at the inlet. As for the other depicted trajectory (blue), the active particle manages to reach the outlet. Along its course through the nozzle it swims upstream several times but in the end the active particle is washed out through the outlet by the background flow. For larger flow rates the trajectories of active particles are less curly, since the flow gets more dominant, see insets of Fig. 5(a).
Next we present results of the second modeling approach which is based on the Yukawasegment model. So far we have concentrated on fixed widths of the inlet and outlet. Here we consider nozzles with fixed length L and inlet width w in and vary nozzle ratio k. We study the behavior of active rods with varied aspect ratio a. As shown in Fig. 7, neither the aspect ratio a, see Fig. 7(a), nor the nozzle ratio k, see Fig. 7(b), have a significant impact on the probability P out which measures how many active rods leave the nozzle at the outlet. However, the aspect ratio a is important for the location where the active rods leave the nozzle at the inlet and the outlet, see Fig. 8. For short rods (a = 2) and small inverse Stokes numbers (σ ≤ 1) the distribution of active particles shows just a single peak located at the center. This peak broadens if the inverse Stokes number increases, which is in perfect agreement with the results obtained by the first approach, cf. Fig. 3. It is more likely for short rods than for long ones to be expelled towards the center due to rotational diffusion. Hence the distribution of particles at the outlet for long rods (a = 10) shows additional peaks close to the wall. These peaks become smaller if the inverse Stokes number increases. The distribution of particles leaving the nozzle at the inlet is similar to our first approach. While the distribution is almost flat for small inverse Stokes numbers, increasing this number makes Figure 7: (a) Probability for an active particle to reach the outlet of the nozzle P out as a function of inverse Stokes number σ for three given aspect ratios a of self-propelled rods and (b) for a fixed aspect ratio a and three given ratios of the nozzle k. Insets show close-ups. it impossible to leave the nozzle close to the center at the inlet. Similar to the outlet the wall accumulation at the inlet is more pronounced for longer rods. By comparing the orientation of the particles at the outlet, the influence of the actual length of the rod becomes visible, see Fig. 9. As seen before for short rods, a = 2, for small inverse Stokes numbers σ there is no wall accumulation. Hence most particles leave the nozzle close to the center and are orientated in the direction of the outlet. This profile smears out if the inverse Stokes number is increased to σ = 1. For larger inverse Stokes numbers the figures are qualitatively similar to the one obtained by the first approach, cf. Fig. 4(c). Particles in the bottom half of the nozzle tend to point upwards and particles in the top half tend to point downwards. The same tendency is seen for long rods a = 10 and small inverse Stokes number. However for long active rods, this is because they slide along the walls. The bright spots close to the walls for long rods and large inverse Stokes numbers indicate that particles close to the walls are flushed through the outlet by the large background flow even if they are oriented upstream. In addition, there are blurred peaks away from the walls for large inverse Stokes numbers σ. The corresponding particles crossed the outlet with mostly upstream orientations. This is similar to Fig. 4(c), where particles exhibiting in-bulk rheotactic characteristics were observed at the outlet of the nozzle.
By comparing the results for individual active rods, see again Fig. 9, with those for interacting active rods at a finite packing fraction ρ = 0.1, see Fig. 10, we find that wall accumulation becomes more pronounced. Mutual collisions of the rods lead to a broader distribution of particles. For long rods, a = 10, the peaks at ϕ ≈ 0 and ϕ ≈ ±π remain close to the walls and the blurred peaks at the center vanish. Figure 10: Outlet distribution histograms for (y, ϕ)| out computed for given inverse Stokes numbers σ and a nozzle with an outlet width w out = 50λ for active rods with an aspect ratio (top row) a = 2 and (bottom row) a = 10 for a packing fraction ρ = 0.1.

Optimization of focusing
Here we study the properties of the active particles in more detail and provide insight into the nozzle geometry, the background flow and the size of the swimmers that should be used in order to optimize the focusing at the outlet of the nozzle.
For this purpose we study three distinct quantities. The averaged dwell time T , the time it takes for an active particle to reach the outlet, the mean alignment of the particles measured by cos ϕ out and the mean deviation from the center y = 0 at the outlet |y out | . As depicted in Fig. 5, for increasing inverse Stokes number the probability for active particles to reach the outlet increases. However they are spread all over the outlet. This is quantified by the |y out | . Small values of |y out | correspond to a better focusing. If particles leave the nozzle with no preferred orientation, their mean orientation vanishes, cos ϕ out = 0; in case of being orientated upstream we obtain cos ϕ out = −1 and finally cos ϕ out = 1 if the particles are pointing in the direction of the outlet. Obviously in an experimental realization a fast focusing process and hence small dwell times T would be preferable.
The numerical results obtained by the first modeling approach are depicted in Fig. 11. While the dwell time hardly depends on the size ratio k of the nozzle, obviously the strength of the background flow has a huge impact on the dwell time and large inverse Stokes numbers σ lead to a faster passing through the nozzle of the active particles, see Fig. 11(a). The alignment of the active particles, cos ϕ out , becomes better if the nozzle ratio k is large and the flow is slow, see Fig. 11(b). The averaged deviation from the centerline |y out | increases with increasing nozzle ratio k since the width of the outlet becomes larger. As could already be seen in Fig. 3, the averaged deviation from the centerline is non-monotonic as a function of the inverse Stokes number and shows the smallest distance from the centerline for all nozzle ratios if the strength of the flow is comparable to the self-propulsion velocity of the swimmers, σ = 1. Let us now study how these three quantities depend on the aspect ratio of the swimmer. To this end, we use the second modeling approach. We consider all three parameters as a function of the inverse Stokes number σ. Longer rods have a shorter dwell time so that they reach the outlet faster, see Fig. 12(a). Increasing the flow velocity obviously leads to a decreasing dwell time. The same holds for the mean alignment -it decreases for increasing inverse Stokes number, see Fig. 12(b). Moreover, for small inverse Stokes numbers, σ ≤ 2, the mean alignment is better for long rods. For large inverse Stokes numbers, long rods a = 10 are washed out with almost random orientation, however short rods a = 2 are slightly aligned with the flow. Short rods are focused better for small inverse Stokes numbers, σ ≤ 2, see Fig. 12(c), due to wall alignment and wall accumulation of longer rods. For larger inverse Stokes numbers, it is the other way around -long rods are better focused. Comparing various nozzle ratios k with fixed simmers' aspect ratio a, we obtain that smaller ratios k lead to smaller dwell times [ Fig. 12(d)] and better alignment [ Fig. 12(e)]. For narrow outlets (small k) the active particles leave the outlet closer to the center, see Fig. 12(f).

Discussion
We discuss the stability of particles around the centerline y = 0 in the presence of a background flow and confining walls if they are converging with a non-zero slope α. This stability is in contrast to a channel with parallel walls, where an active particle swims away from the centerline provided that its orientation angle ϕ is different from nπ, n = 0, ±1, ±2, . . .. Indeed, in the case of a straight channel, α = 0, the background flow is defined as u x = u 0 (H 2 − y 2 ), u y = 0 (Poiseuille flow; u 0 is the strength of the flow, 2H is the distance between the walls). Then the system (3)-(5) reduces tȯ Here we omit the equation for x(t) due to invariance of the infinite channel with respect to x and neglect orientation fluctuations, that is D r = 0. The phase portrait for this system is depicted in Fig. 13(a). Dashed vertical lines ϕ = nπ, n = 0, ±1, ±2, . . . consist of stationary solutions: if an active particle is initially oriented parallel to the walls, it keeps swimming parallel to them. If initially ϕ is different from nπ, then the active particle swims away from the centerline, y(t) → ±∞ as t → ∞. When the walls are converging, α > 0, the y-component of the background flow is nonzero and directed towards the centerline. For the sake of simplicity we take u y = −αy, α > 0 and u x as in the Poiseuille flow, u x = u 0 (H 2 − y 2 ). In this case, the system (3)-(5) reduces tȯ ϕ = −(α/2) sin 2ϕ + u 0 y(1 − cos 2ϕ) (9) y = −αy + v 0 sin ϕ.
The corresponding phase portrait for this system is depicted in Fig. 13(b). Orientations ϕ = nπ represent stationary solutions only if y = 0. In contrast to the Poiseuille flow in a straight channel, see Eqs. (7) and (8), these stationary solutions (ϕ = πn, y = 0) are asymptotically stable with a decay rate α (recall that α is the slope of walls). In addition to these stable stationary points there are pairs of unstable (saddle) points with non-zero y (provided that v 0 > 0). In these saddle points, the distance from centerline |y| does not change, since a particle is oriented away from centerline, so the propulsion force moves the particle away from the centerline and this force is balanced by the convergent component of the background flow, u y , moving the particle toward the centerline. The orientation angle ϕ does not change since the torque from the Poiseuille component of the background flow, u x , is balanced by the torque from the convergent component, u y . We also draw the phase portrait for the converging flow u BG = (u x , u y ) introduced in Section 2, Fig. 13(c). One can compare the phase portraits Fig. 13(b) and Fig. 13(c) around the stationary point (ϕ = 0, y = 0) to see that the qualitative picture is the same: this stationary point is stable and it neighbors with two saddle points.
The asymptotic stability of (ϕ = 0, y = 0) means that if a particle is close to the centerline and its orientation angle is close to 0 (particle is oriented towards the outlet), it will keep swimming at the centerline pointing toward the outlet, whereas in Poiseuille flow the particle would swim away. The asymptotic stability of (ϕ = ±π, y = 0) is evidence of that in the converging flow there is rheotaxis not only at walls but also in the bulk, specifically at the centerline. Another consequence of this stability is the reduction of effective rotational diffusion of an active particle in the region around the centerline, that is the mean square angular displacement ∆ϕ 2 is bounded in time due to the presence of restoring force coming from the converging component of the background flow (cf. diffusion quenching for Janus particles in [48]). Finally, we note that the nozzle has a finite length L and thus, the conclusions of the stability analysis are valid if the stability relaxation time, 1/α s, does not exceed the average dwell time T . We introduce a lower boundT for the dwell time T as the dwell time of an active particle swimming along the centerline oriented forward, ϕ = 0: T = Lk/(σv 0 (1 − k)) ln |1 + σ(1 − k)/(k(σ + 1))|.
Our numerical simulations show thatT underestimates the average dwell time by a factor larger than two. Using this lower bound, we obtain the following sufficient condition for stability: kw in σv 0 ln 1 + σ(1 − k) k(σ + 1) ≥ 1.

Conclusion
In this work we study a dilute suspension of active rods in a viscous fluid extruded through a trapezoid nozzle. Using numerical simulations we examined the probability that a particle leaves the nozzle through the outlet -which is the result of the two counteracting phenomena. On the one hand, swimming downstream together with being focused by the converging flow increases the probability that an active rod leaves the nozzle at the outlet. On the other hand, rheotaxis results in a tendency of active rods to swim upstream. Theoretical approaches introduced in this paper can be used to design experimental setups for the extrusion of active suspensions through a nozzle. The optimal focusing is the result of a compromise. While for large flow rates it is very likely for active rods to leave the nozzle through the outlet very fast, their orientation is rather random and they pass through the outlet close to the walls. The particles are much better aligned with the flow for small flow rates and focused closer to the centerline of the nozzle, however the dwell time of the particles becomes quite large. Based on our findings the focusing is optimal if the velocity of the background flow and the self-propulsion velocity of the active rods are comparable. To reduce wall accumulation, the rods should have a small aspect ratio.
We find that rheotaxis in bulk is possible for simple rigid rodlike active particles. We also established analytically the local stability of active particle trajectories in the vicinity of the centerline. This stability leads to the decrease of the effective rotational diffusion of the active particles in this region as well as the emergence of rheotaxis away from walls. Our findings can be experimentally verified using biological or artificial swimmers in a converging flow.