Dynamical density functional theory for circle swimmers

The majority of studies on self-propelled particles and microswimmers concentrates on objects that do not feature a deterministic bending of their trajectory. However, perfect axial symmetry is hardly found in reality, and shape-asymmetric active microswimmers tend to show a persistent curvature of their trajectories. Consequently, we here present a particle-scale statistical approach of circle-swimmer suspensions in terms of a dynamical density functional theory. It is based on a minimal microswimmer model and, particularly, includes hydrodynamic interactions between the swimmers. After deriving the theory, we numerically investigate a planar example situation of confining the swimmers in a circularly symmetric potential trap. There, we find that increasing curvature of the swimming trajectories can reverse the qualitative effect of active drive. More precisely, with increasing curvature, the swimmers less effectively push outwards against the confinement, but instead form high-density patches in the center of the trap. We conclude that the circular motion of the individual swimmers has a localizing effect, also in the presence of hydrodynamic interactions. Parts of our results could be confirmed experimentally, for instance, using suspensions of L-shaped circle swimmers of different aspect ratio.


Introduction
On the scales of active colloidal particles and self-propelled biological microswimmers [1][2][3][4][5][6][7][8][9][10], thermal fluctuations and other perturbations play a prominent role. They lead to continuous reorientation of the self-propelling objects and therefore to stochastically shaped trajectories [2,3,11]. Even more extreme events are given by stochastic run-andtumble motions. For instance, certain bacteria or alga cells are observed to stop their migration, reorient basically on the spot, and then continue their propulsion [1,12]. Such events lead to kinks on the trajectory. The statistics of both types of buckled motion has been studied in detail, both in experiment and in theory [1,3,11,[13][14][15][16][17][18][19]]. Yet, in the absence of any noise, fluctuations, and perturbations, the self-propelling agents considered in most theoretical analyses would show a deterministic straight motion.
Here, we concentrate on active microswimmers that feature a different behavior. Their individual trajectories are systematically curved. Such a situation can arise only, if for each swimmer the axial symmetry around its propulsion direction is broken.
On the one hand, the symmetry breaking can be induced from outside. For instance, if microswimmers are exposed to local surrounding shear flows, the rotational component of the fluid flow can couple to the orientation of the suspended swimmer [15,[20][21][22][23][24][25]. Continuous reorientation of the propulsion direction leads to curved trajectories. Similarly, the symmetry is broken in the presence of a nearby surface.
If during propulsion a swimmer shows rotations of its body around its axis, these rotations can on one side hydrodynamically interact with the surface. Via such hydrodynamic surface interactions the self-rotation couples to the propulsion direction and the trajectory bends. Also steric interactions can support or induce the effect. Thus circular trajectories are observed for many sperm cells and bacteria close to a substrate [26][27][28][29][30].
On the other hand, real swimmers often bring along a broken axial symmetry by themselves [31]. Hardly any object is really perfectly axially symmetric in shape. On purpose, L-shaped active microswimmers have been fabricated and their persistently curved trajectories were analyzed [32][33][34][35]. If the trajectories, including their persistent bending, are confined to a plane, then circular paths arise. This is what we understand by circle swimmers [36]. Apart from that, for deformable self-propelled particles and self-propelled nematic droplets, the symmetry breaking in shape or structure may also occur spontaneously [37][38][39]. Moreover, imperfections in the self-propulsion mechanism can lead to the symmetry breaking and thus to bent trajectories. An example are cells of the algal Chlamydomonas reinhardtii. If one of the two beating flagella generating self-propulsion is weaker or absent, the cellular paths curve [40,41]. Apart from that, near surfaces bent self-propelled objects tend to follow circular trajectories [42,43]. In modeling approaches, circle swimmers are often realized by simply imposing an effective torque or rotational drive in addition to the self-propulsion mechanism [15,31,35,[44][45][46][47][48][49][50][51][52][53][54][55][56][57][58][59].
We have mentioned above that studies on circle swimmers are relatively rarely encountered when compared to the numbers of works on objects propelling straight ahead. Even less frequent are studies on the collective behavior of circle swimmers [29,43,50,52,54]. Particularly, this applies when hydrodynamic interactions in crowds of suspended microswimmers are to be included.
In equilibrium, i.e., for passive systems, density functional theory [81][82][83][84][85] is, in principle, an exact theory. It can be extended to overdamped relaxational dynamics in terms of DDFT [85][86][87][88] by assuming at each instant an effective equilibrium situation to evaluate the involved potential interactions. For example, solidification processes are addressed in this way [71,72,74,78]. Since microswimmers by construction operate at low Reynolds numbers [89], their dynamics is overdamped. This makes DDFT a promising candidate to study their statistical behavior.
Extending density functional theory to intrinsically non-equilibrium systems, DDFTs for "dry" self-propelling agents had already been developed before [90][91][92] and tested against agent-based simulations [90,92]. Moreover, to characterize passive colloidal particles in suspensions, hydrodynamic interactions had been incorporated into DDFT [93][94][95][96][97][98][99] and agreement was found with explicit particle-based simulations [93,94,96,97]. Our recently developed DDFT for microswimmers incorporates and combines all the central previous ingredients, i.e., self-propulsion, steric interactions between the swimmers, hydrodynamic interactions between the swimmers, as well as exposure to and confinement by external potentials [66]. As we have demonstrated and as further detailed below, this dynamical theory qualitatively reproduces previous simulation results [100,101] in which combined action of all these ingredients determines the overall behavior.
Here, we proceed by an additional step forward. We extend our microscopic statistical characterization (DDFT) to circle swimmers. In this way, we can now characterize the collective behavior of such non-straight-propelling agents, including the effect of hydrodynamic interactions. Only then, for instance and as we will show below, can the symmetry breaking induced by hydrodynamic interactions in a radial confinement be described qualitatively correctly.
We first introduce our minimal model for circle swimmers in section 2. Next, in section 3, we list the extension of the theory. It is evaluated numerically in section 4 to study the behavior of circle swimmers under radial confinement as a function of the bending of their trajectory. A short summary and outlook are given in section 5. a σ/2 Figure 1. Minimal model of a circle microswimmer. Two active force centers separated by a distance L exert the active forces ±f = ±fn onto the fluid, wheren marks the direction of the principal swimmer axis. The resulting fluid flow is indicated by the small light arrows that are rescaled to identical length for visualization. A spherical swimmer body of effective hydrodynamic radius a is, in general, asymmetrically placed between the force centers. Its shift alongn out of the symmetry center, quantified by α = 1/2, leads to forward (or backward) propulsion. An additional transversal shift, quantified by γ, implies circular trajectories for α = 1/2 and γ = 0. This circle swimmer is biaxial, the secondary axis denoted byû. We impose spherically symmetric steric interactions between the swimmers, with an effective interaction radius σ/2.

Minimal model circle swimmer
As outlined above, our goal is to establish a microscopic DDFT of circle-swimmer suspensions. The term 'microscopic' here refers to the length scales of an individual microswimmer. To base our DDFT on such scales, we need to first introduce an explicit minimal model for a microscopic circle swimmer.
Our statistical theory will apply to (semi)dilute suspensions of microswimmers based on their far-field hydrodynamic interactions. Therefore, a minimal model microswimmer is needed that shows the correct leading-order far-field hydrodynamic fluid flows together with a self-consistent description of its self-propulsion. Yet, at the same time, it must be simple enough to still be efficiently included into the statistical description. Figure 1 represents our corresponding swimmer model. The unit vectorn identifies the principal swimmer axis and orientation.
Any active microswimmer exerts forces onto the surrounding fluid. Amongst them are the spatially distributed active forces that initiate self-propulsion. They are generated, for instance, by the rotation of flagella or beating of cilia [12,14,40,41]. In our model, all these active forces are thought to be gathered and concentrated in one spot. In figure 1, this leads to the active point force −f acting on the surrounding fluid. Instantly, due to the nature of the considered low-Reynolds-number motion, see below [89], for a freely suspended microswimmer all these active forces are balanced by frictional forces distributed over its body. We consider all these counteracting frictional forces to be concentrated in another spot, leading to the point force +f in figure 1. These two spots in general do not coincide, depending on the actual swimmer geometry. Here, they are separated by a distance L, see figure 1. However, since no net force nor torque may act on a freely suspended microswimmer, the two forces ±f need to be of same magnitude but oppositely oriented, located and aligned along a common axis. We may thus parameterize them as ±f = ±fn. They act onto the fluid and set it into motion as indicated by the small arrows in the background of figure 1. In analogy to straight-swimming terminology [4], for f > 0, i.e., the depicted case, we call the object a pusher. For f < 0, we term it a puller.
Next, we place a spherical swimmer body of effective hydrodynamic radius a nearby the two force centers. The whole construct is a rigid object, i.e., the force centers and forces have to rigidly translate and rotate together with the sphere, maintaining their mutual distances and orientations. The role of the sphere is purely to realize selfpropulsion of the whole object. Since all forces exerted by the swimmer onto the fluid have already been concentrated into the two force centers above (ignoring all forces that lead to higher-order contributions to the hydrodynamic far-field) the sphere is considered not to exert any remaining force onto the fluid any longer. Its sole role is to be convected by the self-induced fluid flow, leading to the overall self-propulsion. Unless it is positioned into the exact point of symmetry between the two force centers, a net transport of the swimmer results in the induced fluid flow. For a shift of the sphere alonĝ n out of the symmetry plane between the two force centers, the whole object propels into the direction of one of the two forces. This shift is quantified by the parameter α, with α = 1/2 marking the symmetric configuration.
In addition to our swimmer model in Ref. [66], we now consider an extra shift of the spherical swimmer body into a direction perpendicular ton. The parameter γ quantifies this shift, see figure 1, so that the axial symmetry is broken for γ = 0. Consequently, for α = 1/2 and γ = 0, the swimmer in the absence of any fluctuations starts to circle, as quantified below. Moreover, it is now biaxial, with the additional axis marked asû, see figure 1.
Since we consider the hydrodynamic interactions at a far-field level, we need to hinder the microswimmers from coming too close to each other. Therefore, we consider spherically symmetric soft steric interactions between the swimmer bodies of effective radius σ/2 > [(max{α, 1 − α}) 2 + γ 2 ] 1/2 L to maintain an effective distance between them. Altogether, the whole rigid swimmer object in figure 1 is force-and torque-free, as mandatory for a microswimmer suspended in a bulk fluid, see also Appendix A."

Hydrodynamic interactions
We now consider N identical circle microswimmers suspended in the fluid and use indices i = 1, ..., N to label them. As described above, for f = 0, each circle swimmer sets the surrounding fluid into motion due to its active forces exerted by the active force centers. In addition to that, the swimmer bodies may be subjects to forces F i and torques T i . These may, for instance, be stochastic in nature, result from steric interactions between the circle swimmers, or be imposed from outside. Since the dynamics of microswimmers is usually determined by low Reynolds numbers [89], it is described by the linear Stokes equation [102]. That is, their dynamics is overdamped, and the forces F i and torques T i are directly transmitted to the surrounding fluid, setting it into motion. Moreover, since the swimmers are suspended in the fluid, they are translated and rotated by the induced fluid flows. The instantly resulting velocities v i and ω i are calculated from a matrix equation as [66] v for i = 1, ..., N . In (1), the first product on the right-hand side includes the influence of the passive swimmer bodies. µ tt ij , µ tr ij , µ rt ij , and µ rr ij are the familiar mobility matrices that express how swimmer i is translated and rotated due to the forces and torques transmitted by the swimmer body j onto the fluid [66,[102][103][104]. These expressions are the same as for suspended passive colloidal particles and result from an expansion in the inverse separation distance between the swimmer bodies, where here we proceed up to the third order, i.e., the Rotne-Prager level.
The second product on the right-hand side of (1) arises because of the active forces that the swimmers exert onto the fluid. Naturally, these actively induced fluid flows likewise contribute to the velocities v i and angular velocities ω i of all swimmer bodies. Λ tt ij and Λ rt ij are the corresponding mobility matrices. The entries  in these expressions arise because our swimmers do not exert active torques onto the suspending fluid.
More precisely, the mobility matrices Λ tt ij and Λ rt ij describe how the active forces exerted by the two force centers of swimmer j onto the fluid influence the velocity v i and angular velocity ω i of swimmer i, respectively. Since swimmer j carries two active force centers exerting the two forces ±f j = ±fn j , both Λ tt ij and Λ rt ij split into two contributions [66], In contrast to the passive swimmer bodies, the active force centers are point-like. Therefore, the expressions for the four mobility matrices µ tt± ij and µ rt± ij are slightly modified when compared to the corresponding expressions for the hydrodynamic interactions between the passive swimmer bodies in (4) and (6) [66], Here, r ± ij are the distance vectors between the passive body of swimmer i and the active force centers of swimmer j, exerting the forces ±f j = ±fn j onto the fluid, respectively. Again, r ± ij = |r ± ij | andr ± ij = r ± ij /r ± ij . In contrast to Ref. [66], where straight-propelling microswimmers were investigated, we here must take into account the additional transversal shift of the active force centers with respect to the swimmer bodies, see figure 1. Therefore, we now obtain Naturally, the values of α and γ must assure that the force centers of each swimmer are located outside the hydrodynamic radius a of the swimmer body, i.e., [(min{α, 1 − α}) 2 + γ 2 ] 1/2 L > a. We consider our spherical swimmer body to exclusively act as a probe particle. Therefore our active mobility matrices for interactions between different swimmers (i = j) are given to lowest order in (α − 0.5) and/or a/L. In a next step, the distortion of the flow field by the rigid swimmer body could be included by considering the image system within a rigid sphere [105,106].
Moreover, for i = j, (1) together with (7) Radius R s of the circular trajectory at vanishing fluctuations for one single isolated circle microswimmer as introduced in figure 1 (L/a = 3). The color map indicates R s as a function of α and γ. For γ → 0, R s /a → ∞ and the swimmer self-propels straight ahead. We do not allow a force center to be placed within the hydrodynamic radius a from the center of the swimmer body, reflected by the white area on the top right.
speed v s and constant angular speed ω s along a closed circular trajectory of radius R s = v s /ω s for all times. Since both v s and ω s depend on the relative position between the swimmer body and the two force centers, R s can smoothly be tuned between almost zero and infinity by altering the parameters α and γ; see figure 2. Moreover, both Λ tt ii and Λ rt ii are independent of f . Thus, both v s and ω s scale linearly with f , see (7)- (12). Therefore, R s is independent of the active force f . Swimming faster does not change the radius of the circle.
Technically, our mobility matrices represent the solutions to the underlying Stokes equation for the flow of the suspending fluid at low Reynolds number [102]. In this way, the role of the fluid is implicitly included in our description.

Stochastic forces, external forces, and steric interactions
Our remaining task is to specify the forces F i and torques T i acting on the swimmer bodies in (1). The forces are set to Here, the first contribution represents the effective influence of the stochastic forces due to thermal fluctuations [107]. k B is the Boltzmann constant, T the temperature, is the probability density to find at a certain time t the swimmers at positions r i with orientationsn i andû i , i = 1, . . . , N . From this form, the correct diffusional behavior is reproduced in the statistical approach, see below.
The overall potential in the second part of (13) reads In this expression, the first term describes the steric interactions between the swimmer bodies. We here choose a soft GEM-4 potential of the form [78,108] where 0 sets the strength of the interactions. u ext is an external potential acting on each swimmer body and further addressed below. Finally, the only torques that we consider to act on our spherically symmetric swimmer bodies are stochastic ones, Here, the operator ∇ or i contains the derivatives with respect to the particle orientations. If the swimmers and their orientations are confined to a flat plane, for instance, the xy plane in Cartesian coordinates, one angle ϕ i is sufficient to characterize the orientation of each swimmer i. Then the operator reduces to ∇ or i =ẑ ∂/∂ϕ i , whereẑ is the (oriented) Cartesian unit vector perpendicular to the xy plane in a three-dimensional Euclidean space. In three dimensions, explicit expressions using Eulerian angles exist [91,109].

Dynamical density functional theory for circle swimmers
Based on our minimal microswimmer model, we can now derive a microscopic statistical description in terms of a dynamical density functional theory (DDFT) for suspensions of identical circle swimmers. The derivation follows the same lines as in our previous work on straight-propelling microswimmers [66]. However, several changes result from the present biaxiality of the individual swimmers.
We start from the microscopic Smoluchowski equation which states the conservation of the overall probability density. Here, we have to insert the swimmer velocities v i and angular velocities ω i as given by (1)- (16). Although v i and ω i depend on ln P via (13) and (16), it is important to stress that (17) is still linear in P . Using the chain rule in (13) leads to ∇ i ln P = (∇ i P )/P , which in combination with the factor P in (17) leads to the linear contribution ∇ i P . The same argument applies to the term ∇ or i in (16) when inserted into (17). To obtain from (17) the n-swimmer density of finding n of the identical N circle swimmers at a certain time at certain positions with certain orientations, we must integrate out from (17) all but the degrees of freedom of n swimmers. We denote by X i all degrees of freedom of the ith swimmer. Then, the n-swimmer density is obtained from the overall probability density P as In the special case of all swimmers and their orientations being confined to a flat plane, Our goal is to obtain an equation for the dynamics of the one-swimmer density ρ (1) (X, t) to find a circle swimmer at time t with position and orientation(s) X. However, the integration scheme in (18) leads to a non-closed equation for the time derivative of ρ (1) . Because of our pairwise hydrodynamic and steric interactions, ρ (1) couples to the pair density ρ (2) , and, in combination of both interactions, also to ρ (3) [66,93,94]. This starts a whole hierarchy of coupled dynamical equations, called BBGKY hierarchy [81]. To close the dynamical equation for ρ (1) , we need to express the densities ρ (2) and ρ (3) in this equation as a function of ρ (1) . DDFT provides a strategy by mapping each state of the system instantaneously to a corresponding equilibrium situation [85][86][87][88].
For this purpose, we recall that an external potential enters the dynamical equation via (14). At each moment in time, DDFT assumes that the instant state of the system is caused by an effective external potential Φ ext . This Φ ext intermittently takes the place of our physical external potential u ext .
In equilibrium, density functional theory (DFT) implies that Φ ext is uniquely determined by the density ρ (1) [81][82][83][84][85][86][87][88]. It follows by minimizing the grand canonical potential functional Ω with respect to ρ (1) . Here, is the entropic free-energy functional for ideal non-interacting particles, with λ the thermal de Broglie wave length [110]. An exact expression for the excess free-energy functional F exc [ρ (1) ], which contains all particle interactions beyond the idealized noninteracting limit, is typically not known and needs to be approximated. The third functional describes the interactions with the external potential, where the effect of a chemical potential is implicitly included into Φ ext . Minimizing Ω with respect to ρ (1) leads to the equilibrium relation In equilibrium the swimmers are inactive (f = 0). Then, we may further argue that the corresponding N -swimmer probability density P eq solely depends on the overall potential U = U (X 1 , . . . , X N ) as in (14), but with Φ ext taking the place of u ext . Then, P eq should follow the Boltzmann form with β = (k B T ) −1 . Applying to this relation the positional gradient for the ith swimmer, we obtain We then follow (18) and integrate out all coordinates from this relation except for those of the ith swimmer. Since all swimmers are identical, this leads to the so-called YGB relations of first order [81,109], The YGB relations of second order are obtained by integrating out from (24) all coordinates but those of the ith and one other swimmer [81,109], resulting in We then eliminate Φ ext from the last two equations by inserting (22). The resulting relations dX ρ (2) (X, X )∇ r u(r, r ) = ρ (1) (X)∇ r δF exc δρ (1) (X) (27) and have the same structure as the corresponding ones in Ref. [66]. DDFT assumes that these relations are still instantly satisfied in non-equilibrium at each moment in time. All contained quantities are then assumed to be dynamical and non-equilibrium ones. In this way, they are inserted into the dynamical equation for ρ (1) . Our assumption implies that the higher-order swimmer densities relax quickly when compared to the lower-order ones [111]. Since our motion at low Reynolds numbers is overdamped, it is conceivable that this adiabatic approximation leads to reasonable results. Previous comparison with particle simulations has confirmed this assertion qualitatively [66].
Altogether, we obtain from this procedure where J 1 , . . . , J 6 are current densities.
They are of similar structure as the corresponding quantities in Ref. [66], but particularly the active current densities J 3 and J 6 differ in the present case because of the transversal shift of the active force centers, see figure 1, We note that, in our case, J 2 = 0 in (31), and also the integral containing ρ (2) in (34) vanishes. The reason is the spherical shape of our passive swimmer bodies, resulting in hydrodynamic interactions that do not depend on the swimmer orientations.
Here, V (r, r ) = u(r, r ), if r = r . For r − r → , we let βV → ∞ to avoid the hydrodynamic divergence that appears in the unphysical situation of two swimmers being located at the same position. This corresponds to exp(−β 0 ) → 0, which in our typical choice of parameters will be exp(−10) already without setting exp(−βV ) to zero at r = r . In this way, our dynamical equation for ρ (1) is derived and finally closed.
To demonstrate the power of our DDFT for circle swimmers, we now address the confinement in a spherically symmetric trap. In particular, we focus on the effect of an increasing curvature of the swimming paths.

Circle swimmers in a spherically symmetric trap
To address planar geometries, we confine the center of mass of each swimmer i as well as its two orientation vectorsn i andû i to the Cartesian xy plane so thatn i ×û i =ẑ. Then, one angle ϕ i is sufficient to parameterize the swimmer orientation, see our remarks below (16). We measure ϕ i relatively to the x axis. Still, three-dimensional hydrodynamic interactions apply. One possible realization of this geometry are swimmers confined to the interface between two immiscible fluids of identical viscosity.
Next, we specify the spherically symmetric confining external potential in (14). As in Ref. [66], we use a quartic potential centered in the origin, where r k = |r k |. This potential is more shallow around the center and then shows a steeper increase than a harmonic trap, which partially emphasizes the effects that we address in the following. Yet the precise functional form is not relevant for their qualitative nature. To evaluate our DDFT numerically, the finite-volume-method (FVM) partialdifferential-equation solver FiPy [113] is employed. Our numerical grid is regular, quadratic in the xy space, and typically consists of 80 × 80 × 16 grid points for the x, y, and ϕ coordinates, respectively. (Non-orthogonal meshes might produce significant numerical errors due to the assumption of orthogonality by the solver [114,115]. We avoid this by using an orthogonal grid.) We only analyze the behavior in one single isolated trap. Nevertheless, numeric periodic boundary conditions are imposed in all directions for technical reasons to allow for Fast Fourier Transformation. To avoid unphysical feedback between particles through the walls of the box, long-ranged hydrodynamic interactions are cut at distances larger than half a box length. Care is taken that the extension of the density cloud, before it basically decays to zero due to the external potential, is smaller than half a box length. In this way, the density cloud does not interact with itself through the periodic box boundaries. However, the box is large enough to account for all hydrodynamic interactions within the effective confinement by the spherical trap. The steric interactions in (15) do not need to be cut as they quickly decay with increasing distance. If, instead of one single isolated trap, an array of periodically placed traps were to be regarded, one would have to account for the long-ranged hydrodynamic interactions between the individual traps including the periodic images of the system, e.g., via Ewald summation techniques [116][117][118].  To display our results, we extract the spatial swimmer density ρ(r, t) = dϕ ρ (1) (r, ϕ, t) (39) and the orientational vector field n (r, t) = dϕn(ϕ) ρ (1) (r, ϕ, t) from our calculations. These two fields are indicated by color plots and by white arrows, respectively, in the figures referred to below. In these plots, the spatial density ρ(r, t) is normalized by the densityρ averaged over the whole simulation box.
As an initial condition, we start from randomized density distributions. The system is then equilibrated in the trap with self-propulsion switched off, f = 0. The density quickly relaxes into a radially decaying distribution with a small central dip stemming from steric repulsion, see figure 3. We measure time t in units of σ 2 /(µ t k B T ). At t = 0, self-propulsion is switched on. Such a process could be achieved in reality, for instance, using light-activated synthetic swimmers [11,16,34,[119][120][121]. If, for example, activation of self-propulsion is sensitive to the wavelength of the irradiated light [119], confinement might be achieved simultaneously by optical trapping using light of a different frequency.
To characterize the relative strength of self-propulsion, often the dimensionless Péclet number Pe tr is introduced [101]. In our context, it measures the ratio of active to diffusive passive motion on a relevant length scale, here set by σ. Therefore, In the following, we concentrate on microswimmers of significant activity, Pe 1. Moreover, we may in the case of circle swimming analogously define a rotational Péclet number, For Pe rot ≈ 0, the curvature of the swimmer trajectory is negligible. In our numerical scheme, we directly set the parameters determining the geometry of the swimmers in figure 1. The corresponding Péclet numbers can then be extracted by calculating v s and ω s from (1) and (7)- (12) for i = j. It turns out that increasing the character of circle swimming, i.e., decreasing the radius of the unperturbed swimming path, see figure 2, has a qualitative effect on the appearance of the trapped swimmer suspension. To demonstrate this, we first further analyze some results of straight swimming [66] obtained by our modified simulation scheme and then compare with the results for circle swimming.

Straight swimming
Straight motion of the individual swimmers is enforced in our approach by setting γ = 0, see figure 1 [66]. For straight propelling objects under spherical confinement, the formation of a high-density ring has been reported several times [66,100,101,122,123]. In agreement with previous studies, the formation of a high-density ring can be reproduced after switching on the active drive in our simulations. This ring is particularly symmetric when we switch off the hydrodynamic interactions between the swimmers, see figure 4 (a). Its approximate radius is determined by balancing the active forward drive of the swimmers with the confining external potential force, leading to R ring ∼ σ(v s σ/4µ t V 0 ) 1/3 = σ(Pe tr k B T /4V 0 ) 1/3 .
In the next two rows, figure 4 shows the behavior when hydrodynamic interactions between the swimmers are included as prescribed by our DDFT. They have a qualitative impact. The high-density ring at the investigated propulsion strengths develops a tangential instability and the circular symmetry is broken. Also this effect has been described before [66,100,101]. The swimmers tend to polarly order in the emerging high-density spot while propagating against the confining potential. Consequently, they collectively pump the surrounding fluid into the opposite direction. Thus the configuration was referred to as a "hydrodynamic fluid pump" [66,101]. Here, we observe that the effect is stronger in figure 4 (b), which depicts the result for pushers, f > 0. In contrast to that, figure 4 (c) was obtained with the sign of the active forces flipped to f < 0, describing a suspension of pullers, yet with all other parameters unchanged. Obviously, the tangential symmetry breaking is restricted in the latter case.
The cause of this spontaneous symmetry breaking was attributed in Ref.
[101] to a positive feedback mechanism. If a spot of higher density appears on the ring, with the swimmers collectively pushing against the external potential, the resulting oppositely oriented fluid flow rotates nearby swimmers towards the high-density area. Consequently, they join the spot of higher concentration. In our DDFT, this effect is included by the contribution ∼ µ rt r,r ∇ r u ext (r ) to the current density J 4 in (33).  Additionally, pushers actively generate inward flows from their sides, see figure 1. When the swimmers are pointing outward on the ring, this further supports their lateral concentration, see the illustration in figure 5 (a). Here, these active contributions are represented by the second term in the current density J 3 in (32). In contrast to that, for pullers, the actively induced flow fields are inverted. This in effect repels outward pointing swimmers on the ring from each other, see also our schematic illustration in figure 5 (b). The qualitative schematics in figure 5 (c)-(f) indicate that also the curvature of the high-density ring may have a significant impact via the current density J 6 and lead to differences between pushers and pullers. The relative magnitudes of all these different contributions basically involve all system parameters, i.e., temperature T , the viscosity η of the fluid, the strength of the active force f , the nature of the swimmer (pusher vs. puller), the strength and radius of the steric interactions, and the strength of the confinement. After the formation of the high-density spot, see figure 4 (b), at strong enough active drive f > 0, we observe yet another spontaneous symmetry breaking. In the rightmost snapshot of figure 4 (b), the averaged self-propulsion directions do not point radially outward any more. Instead, they have tilted to one side towards the tangent of the previous high-density ring. For straight swimming objects, the selection of one of the two tilting directions depends solely on small variations in the initialization of the system.
As a result of the tilting, the high-density spot starts to circle around the trap, smearing out the faded ring to some extent. Depending on the parameters, we may nearly recover a high-density ring, however, now with the swimmer orientations not pointing outward. Interestingly, for suspensions of pullers at elevated |f |, we so far have not observed this behavior. Instead, again a ring of radially oriented swimmers emerges, see figure 4 (c). It appears approximately in the same way as for the case without hydrodynamic interactions in figure 4 (b). This behavior is in line with our interpretation of the role of the current density J 3 ∼ f of repelling pullers from each other.
We note that the active current density J 6 in (35) has the potential to drive the secondary spontaneous symmetry breaking observed in the rightmost snapshot of figure 4 (b). Comparing the strength of J 6 with the one of J 4 may also explain the initial formation of the high-density spot as a first instability and then the observed secondary instability. First, on the high-density ring, the swimmers on average feature a larger mutual separation, see the second snapshot of figure 4 (b). Then, at these larger interswimmer distances r ij , the contribution in the current density J 4 driving the spot formation scales as ∼ |µ rt ij | ∼ r −2 ij . In contrast to that, in the active current density J 6 , we find a scaling ∼ r −3 ij at large interswimmer distances (the two oppositely oriented active forces of each swimmer together appear as a force dipole at larger distances, which reduces the exponent in the scaling of Λ rt ij to ∼ r −3 ij ). Therefore J 4 dominates and can drive the spot formation. At reduced separation in the high-density spot, the active forces are resolved individually and the influence of J 6 can become substantial Figure 5. Effects that can contribute to the observed varying tendency of forming a high-density spot for pushers and pullers in figure 4 (b) and (c), respectively. Large straight arrows of the same color as the swimmer bodies represent the active forces, while the smaller curved arrows indicate the corresponding flow fields. (a) Due to their actively induced transversal inflow of fluid, see figure 1, pushers laterally attract each other hydrodynamically. This supports spot formation in figure 4 (b) when the swimmers are aligned next to each other on a high-density ring. (b) In contrast to that, pullers laterally repel each other hydrodynamically, which counteracts a spot formation, see figure 4 (c). Both effects are described by the contribution ∼ Λ tt r,r to the current density J 3 in (32). (c)-(f) Different effects are possible for the active rotation-translation coupling described by the contribution ∼ Λ rt r,r to the current density J 6 in (35). (c) For pushers next to each other on a high-density ring of high curvature, the inward pointing force center of one swimmer is in close vicinity of the body of the other swimmer, and vice versa. This leads to actively induced rotations of the swimmers and their propulsion directions towards each other, supporting the formation of a high-density spot. (d) The situation is reversed for pullers, leading to effective rotations away from each other. (e) In contrast to the configuration in (c), for low curvature of the high-density ring, the outward pointing active force of one swimmer is closer to the body of the other swimmer, and vice versa. In this way, the swimmers tend to turn away from each other. (f) Along the same lines, pullers also for low curvature of the high-density ring turn away from each other, again counteracting the formation of a high-density spot.  when comparing with J 4 . Now both scale as ∼ r −2 ij , but for elevated |f | the importance of J 6 grows.

Circle swimming
We now turn to increasingly biaxial swimmers for γ = 0, see figure 1. Depending on the values of both parameters α and γ, the unperturbed individual swimmers then show circular trajectories, see figure 2. Particularly, we analyze the changes in the behavior of the suspension when we stepwise increment γ. For each value of γ, we again start from an equilibrated passive initial system and then switch on self-propulsion at t = 0, as before.
By and large, we do not observe abrupt modifications in the overall behavior. Instead it changes rather gradually with increasing γ. For small γ = 0, the behavior of the straight swimming objects is reproduced qualitatively. Only for pushers of stronger active drive f > 0, we note an illustrative alteration. While the sense of circling of the high-density spot around the trap as illustrated in the rightmost snapshot of figure 4 (b) resulted from spontaneous symmetry breaking for γ = 0 and could be clockwise or counterclockwise, it is now increasingly dictated by the sense of the circular swimming trajectory. A comparison between straight swimmers and weak circle swimmers is included by figure 6.
Remarkably, the overall appearance of the suspension changes qualitatively when the nature of circle swimming becomes more pronounced. In our set of parameters we achieve this by increasing γ. The bending of the swimmer trajectories has a localizing effect, as illustrated in figure 7. There, all snapshots show the long-term behavior of the corresponding suspension. From left to right in each row, the strength of circle swimming grows. Due to their persistent self-rotation, the outward propagation of the swimmers against the confining trapping potential is restricted. As a consequence, the concentration of the swimmers in the center of the trap increases. At high enough γ, the density is again peaked around the center of the trap. Comparing figure 7 (a), where hydrodynamic interactions have been switched off, to figure 7 (b) and (c), we infer that hydrodynamic interactions significantly delay the localization around the center of the trap with increasing γ. Yet, at high enough values of γ (rightmost column in figure 7) the localization dominates in all cases. Comparing pushers and pullers in figure 7 (b) and (c), respectively, we note the more persistent nature of the high-density ring in the case of pullers at smaller values of γ, before the collapse towards the center of the trap occurs.
To quantify the modified appearance of the suspension with increasing γ, we introduce the following order parameters. First, we evaluate where in this expression spatial positions r are parameterized by polar coordinates r = (r, ϑ). K(t) becomes nonzero when a tangential instability occurs that breaks the circular symmetry of a high-density ring, leading to an off-center high-density spot. Next, we define withv s for each swimmer denoting the hypothetical instantaneous unperturbed direction of self-propulsion. For γ = 0,v s points along ±n according to the sign of f , but it becomes slightly tilted towardsû for γ = 0. M r (t) quantifies the overall degree of swimmer orientations along the radial direction.
In analogy to that, to quantify the ordering of the swimmer orientations along one of the two tangential directions, the order parameter is evaluated. In the absence of any local orientational order, both M r (t) and M t (t) vanish. For steady-state systems, all three of the above order parameters no longer depend on time in the long-term limit. Figure 8 shows the long-term values of the order parameters K, M r , and M t with increasing biaxiality and degree of circle swimming γ. When hydrodynamic interactions  For pushers (f > 0), we observe with increasing γ that first the high-density spot smears out to a high-density ring that broadens and for high γ collapses towards the center. (c) For pullers (f < 0) the high-density ring appears a bit more stable at lower values of γ, but again a concentration around the center occurs at high γ.
are switched off, for γ = 0 a high-density ring is formed with the swimmers radially aligned, see figure 4 (a). Therefore, K and M t are low, while M r is high.
Including hydrodynamic interactions, pullers (f < 0) here behave in a very similar way, see also figure 7 (c). In contrast to that, pushers (f > 0) show a concentration in high-density spots for γ = 0, see figures 6 and 7 (b), leading to an elevated value of K. Moreover, the self-propulsion directions in this high-density spot by spontaneous symmetry breaking can lean towards one of the two tangential directions, see figures 4 (b) and 6. Therefore, M r and M t are reduced and elevated, respectively, when compared to the other systems in figure 8.
As the degree of circle swimming increases with γ and the swimmers tilt away Order parameters K, measuring the degree of off-center density concentration in a high-density spot, as well as M r and M t , measuring the degrees of swimmer orientations along the radial and one of the tangential directions, respectively, for the systems in figure 7 with increasing biaxiality parameter γ. Again, the situation without hydrodynamic interactions between the swimmers ("no H.I."), the case of f > 0 ("pusher"), and the case of f < 0 ("puller") are depicted. Generally, with increasing γ, off-center concentration in non-rotationally symmetric structures diminishes (drop of K) and the swimmers tilt away from the radial direction (decreasing M r ). We find smooth intermediate transitions for M r and M t around the value of γ that leads to R s = R ring , as indicated by the vertical gray lines. For the parameter chosen here, K drops to zero at very low biaxiality, which is not true for all parameters, see, e.g., figure 6. from the radial outward direction, M r generally decreases. M t first increases as the orientational order shifts from radial to tangential. It then saturates and again slightly decays for high γ, i.e., for small swimming radii. The latter slow decay is supported by the increasing localization in the center of the trap where orientational order vanishes by the overall rotational symmetry. The smooth changes of M r and M t in figure 8 indicate that the transition from off-center high-density rings or spots to centrally localized distributions with increasing γ is rather continuous. This transition should occur when the radius R s of the unperturbed swimmer trajectories and the characteristic radius of the trap R ring become approximately identical. We have indicated the corresponding value of γ in figure 8 by the vertical gray lines.
To also quantify the depletion of the swimmer density in the center of the trap when high-density rings or off-center high-density spots occur, in contrast to the central accumulation when the localizing effect of circle swimming becomes strong, we introduce additional order parameters for ν = 0 and 1. Here, J ν are the Bessel functions of first kind. By construction, O 0 (t) is large when the density is concentrated in the center of the trap, is elevated for off-center distributions.
As demonstrated by figure 9, the transitions as a function of the biaxiality parameter γ are again smooth. Yet, the increasing localization in the center of the trap for increasing γ is obvious. Particularly in the transitional regime, that is, for intermediate values of γ, hydrodynamic interactions apparently counteract localization in the center of the trap. Moreover, for our set of parameters and at large γ, the central concentration of pushers is slightly lower than the one for pullers, if only the sign of f is inverted and all other parameters are kept the same.

Conclusions
In summary, we have presented a microscopic statistical approach in the framework of dynamical density functional theory (DDFT) on active circle swimmers. Hardly any real microswimmer is a perfectly symmetric straight swimmer. Therefore, investigations on the effect of bent migration trajectories are mandatory.
Our theory captures self-propulsion along swimming paths of different preferred curvature, steric and hydrodynamic interactions between the microswimmers, as well as confinement by an external potential. In contrast to many previous descriptions, the curved motion in our case is not directly imposed by an effective torque or angular frequency on the swimmer body. Here, it naturally follows from the geometric structure of our microscopic minimal swimmer model and resulting hydrodynamic effects.
Persistently bent swimming trajectories reduce the global mobility of the swimmers. To study this localizing effect, we analyzed the behavior of microswimmer suspensions in a circularly symmetric trapping potential for increasing degree of circle swimming. no H.I. pusher puller Figure 9. Same as in figure 8, but for the order parameters O 0 and O 1 that measure the degree of near-center and off-center concentration, respectively. At smaller values of the biaxiality parameter γ, O 0 is low and O 1 is high when a high-density ring or a high-density off-center spot has formed. In contrast to that, high O 0 and low O 1 signal a localization of the density around the center of the trap for pronounced circle swimming at high values of γ (see also the rightmost column in figure 7). Again, the vertical gray lines indicate the value of γ implying R s = R ring . Apparently, hydrodynamic interactions slightly counteract the concentration around the center. For strong circle swimming, pullers appear more concentrated around the center than pushers.
Moreover, we distinguished between pusher and puller circle swimmers, and also studied the effect of hydrodynamics by comparison with switched-off hydrodynamic interactions between the swimmers. Straight swimming objects tend to spread out towards the confinement until their active drive is balanced by the confining potential [66,100,101,122,123]. This leads to high-density rings. Such rings may get unstable due to hydrodynamic interactions, particularly for pusher swimmers, leading to the formation of off-center high-density spots [66,100,101]. We have further investigated and quantified this scenario.
Circle swimming can qualitatively affect the behavior. Increasing the degree of circular self-propulsion supports a persistent circling motion of the high-density spots around the trap. At high degrees of circle swimming, the swimmers become localized around the center of the trap, while hydrodynamic interactions seem to slightly counteract this effective confinement. The transition from the off-center towards the centered density distributions appears to be smooth, and we quantified it by introducing several corresponding order parameters.
A long-term goal to extend the present theory would be the characterization of motility-induced phase separation into a dense clustered state and a surrounding low-density gas-like state [119,120,. This phenomenon was observed in particle-based simulations of active Brownian particles [120, 125, 127, 130, 133-136, 138, 142, 144, 145] and described by different statistical or continuum approaches [126,128,129,131,132,136,137,141]. So far, the effect of hydrodynamic interactions on this scenario has only rarely been addressed [134,138]. Our DDFT by construction contains self-propulsion driving the phase separation, steric interactions to avoid a collapse of the clustered state, and hydrodynamic interactions. In previous theoretical approaches, input for the density dependence of the swimming speed [128] or for the front-back imbalance of the pair-correlation function [126,131,141] was required to capture the phenomenon. An interesting question for statistical theories and DDFT is whether such an input will further be necessary in the future, or whether the theories will provide it in a self-consistent way, as encouraged by a recent theoretical study [146]. Moreover, one could then analyze how the clustering behavior is influenced by the circular swimming paths. Apart from that, in the future also the dynamic behavior of pure active swimming rotors [147][148][149] could be considered in an analogous statistical approach, including the induced hydrodynamic interactions between the rotors. We note that, in a different context, for reorienting the swimming motion, e.g., by external fields the consequences have been analyzed for the translational behavior and for the swim stress and pressure [150,151]. Possibly, the latter quantities could also be extracted using our approach and explicit swimmer model. Another extension concerns the treatment of crystallization effects [78] for active microswimmers including hydrodynamic interactions.