A particle-field representation unifies paradigms in active matter

TL;DR

This paper introduces a unified continuum field approach to active matter, revealing how particle shape and rigidity influence emergent order, bridging different paradigms like phase separation and collective motion.

Contribution

It develops a novel particle-field representation inspired by wave-particle duality, unifying active matter behaviors and demonstrating how stresses and shape control emergent states.

Findings

Anisotropic rods form polar, moving clusters.

Soft rods exhibit nematic, large-scale structures.

Symmetry of order is an emergent property, not dictated by interactions.

Abstract

Active matter research focuses on the emergent behavior among interacting self-propelled particles. Unification of seemingly disconnected paradigms -- active phase-separation of repulsive discs and collective motion of self-propelled rods -- is a major challenge in contemporary active matter. Inspired by the quanto-mechanical wave-particle duality, we develop an approach based on the representation of active particles by smoothed continuum fields. On the basis of the collision kinetics, we demonstrate analytically and numerically how nonequilibrium stresses acting among self-driven, anisotropic objects hinder the formation of phase-separated states as observed for self-propelled discs and facilitate the emergence of orientational order. Besides particle shape, the rigidity of self-propelled objects controlling the symmetry of emergent ordered states is as a crucial parameter: impenetrable, anisotropic rods are found to form polar, moving clusters, whereas large-scale nematic structures emerge for soft rods, notably separated by a bistable coexistence regime. These results indicate that the symmetry of the ordered state is not dictated by the symmetry of the interaction potential but is rather a dynamical, emergent property of active systems. This unifying theoretical framework can represent a variety of active systems: living cell tissues, bacterial colonies, cytoskeletal extracts as well as shaken granular media.

Equations

$${\psi }_k\left(\mathbf{r}\right)=e^{-\frac{{\left\{(\mathbf{r}-{\mathbf{r}}_k)\cdot {\mathbf{e}}_{\text{vrule}depth=0.0pt\text{vrule}depth=0.0pt}\left[{\phi }_k\right]\right\}}^2}{2l_{\text{vrule}depth=0.0pt\text{vrule}depth=0.0pt}^2}-\frac{{\left\{(\mathbf{r}-{\mathbf{r}}_k)\cdot {\mathbf{e}}_{\perp }\left[{\phi }_k\right]\right\}}^2}{2l_{\perp }^2}}.$$

$$\epsilon =\frac{l_{\text{vrule}depth=0.0pt\text{vrule}depth=0.0pt}^2-l_{\perp }^2}{l_{\text{vrule}depth=0.0pt\text{vrule}depth=0.0pt}^2+l_{\perp }^2}$$

$${\mathcal{I}}_{kj}\propto e^{-\frac{\left({\mathbf{r}}_k-{\mathbf{r}}_j\right)\cdot \left[\text{mathds}1-\frac{\epsilon }{2}\left(\mathcal{Q}\left[{\phi }_k\right]+\mathcal{Q}\left[{\phi }_j\right]\right)\right]\cdot \left({\mathbf{r}}_k-{\mathbf{r}}_j\right)}{2\left[1-{\epsilon }^2{\cos }^2\left({\phi }_k-{\phi }_j\right)\right]\left(l_{\text{vrule}depth=0.0pt\text{vrule}depth=0.0pt}^2+l_{\perp }^2\right)}},$$

$$u_2\left(\Delta \mathbf{r};\phi ,{\phi }^{\prime }\right)=\kappa \mathcal{F}\left[e^{-\frac{\Delta \mathbf{r}\cdot \left[\text{mathds}1-\frac{\epsilon }{2}\left(\mathcal{Q}\left[\phi \right]+\mathcal{Q}\left[{\phi }^{\prime }\right]\right)\right]\cdot \Delta \mathbf{r}}{2\left[1-{\epsilon }^2{\cos }^2\left(\phi -{\phi }^{\prime }\right)\right]\left(l_{\text{vrule}depth=0.0pt\text{vrule}depth=0.0pt}^2+l_{\perp }^2\right)}}\right],$$

$$m_2\left(\mathbf{r},\phi ,{\phi }^{\prime }\right)=$$

$$+$$

$${\dot{\mathbf{r}}}_k$$

$${\dot{\phi }}_k$$

$${\mathbf{\eta }}_k(t)=\sqrt{2D_{\text{vrule}depth=0.0pt\text{vrule}depth=0.0pt}}{\mathbf{e}}_{\text{vrule}depth=0.0pt\text{vrule}depth=0.0pt}\left[{\phi }_k\right]{\eta }_{\text{vrule}depth=0.0pt\text{vrule}depth=0.0pt,k}(t)+\sqrt{2D_{\perp }}{\mathbf{e}}_{\perp }\left[{\phi }_k\right]{\eta }_{\perp ,k}(t).$$

$${\partial }_{t}P=$$

$$\displaystyle-\!\partial_{\varphi}\big{[}\mu_{\varphi}MP\Big{]}+D_{\varphi}\partial_{\varphi}P,$$

$$\mathbf{F}$$

$$M$$

$$\displaystyle m_{j\rightarrow k}\propto\varepsilon\sin\!\Big{[}2\big{(}\varphi_{k}\!\!\>-\arg\!\left[\mathbf{r}_{k}-\mathbf{r}_{j}\right]\!\!\>\big{)}\!\!\>\Big{]},$$

$$\mathbf{F}$$

$$M$$

$${\partial }_t\rho$$

$${\partial }_t\mathbf{p}$$

$$\delta \mathbf{p}\simeq -\frac{(v_0-2\bar{\mu }{\zeta }_0\kappa {\rho }_0)+\epsilon {\mu }_{\phi }{\zeta }_1\kappa {\rho }_0}{2D_{\phi }}\nabla \delta \rho .$$

$$\Gamma =D_0+\frac{\left(v_0-\bar{\mu }{\zeta }_0\kappa {\rho }_0\right)\left[\left(v_0-2\bar{\mu }{\zeta }_0\kappa {\rho }_0\right)+\epsilon {\mu }_{\phi }{\zeta }_1\kappa {\rho }_0\right]}{2D_{\phi }}.$$

$$\epsilon >\frac{v_0-4\sqrt{D_0{D}_{\phi }}}{{\mu }_{\phi }{\rho }_0{\zeta }_1\kappa }.$$

$$\dot{\mathbf{p}}$$

$$\dot{\mathfrak{Q}}$$

$${\sigma }_p$$

$${\sigma }_n$$

$$g_2^{(0)}$$

$$=\frac{1}{2\pi }{\int }_{-\pi }^{\pi }d\phi g_2\left(\left|\mathbf{r}-{\mathbf{r}}^{\prime }\right|,\arg \left({\mathbf{r}}^{\prime }-\mathbf{r}\right),\phi \right)$$

$$g_2^{(1c)}$$

$$=\frac{1}{\pi }{\int }_{-\pi }^{\pi }d\phi \cos \phi g_2\left(\left|\mathbf{r}-{\mathbf{r}}^{\prime }\right|,\arg \left({\mathbf{r}}^{\prime }-\mathbf{r}\right),\phi \right)$$

Full text

I Introduction

Interacting self-propelled particles, termed active matter, are the standard model of collective behavior out of thermodynamic equilibrium Romanczuk et al. (2012); Marchetti et al. (2013); Menzel (2015). Active systems become increasingly popular in different disciplines studying diverse systems from micron-sized synthetic swimmers Paxton et al. (2004); Palacci et al. (2013); Theurkauff et al. (2012) or phoretic colloids Buttinoni et al. (2013); Golestanian et al. (2007), collective swimming and self-organization of bacteria Dombrowski et al. (2004); Sokolov et al. (2007) as well as self-assembly in biomimetic systems Sanchez et al. (2012) to the behavior of insect swarms Attanasi et al. (2014), sheep herds Ginelli et al. (2015), fish schools Ward et al. (2008) and flocks of birds Ballerini et al. (2008).

The current theoretical understanding of active matter is based on two cornerstones. One of these paradigms is the emergence of phase-separated states in ensembles of self-propelled discs Fily and Marchetti (2012); Redner et al. (2013); Digregorio et al. (2018), also referred to as active Brownian particles 111The term active Brownian particles, originally introduced in Schimansky-Geier et al. (1995) referring to Brownian particles with the ability to generate a field which influences their motion, is often used synonymously in a more general context with self-propelled motion far from equilibrium Schweitzer et al. (1998); Romanczuk et al. (2012)., due to the combined effect of self-propulsion and isotropic repulsion Cates and Tailleur (2015). This phenomenon, sharing similarities with granular media agitated by vibration Aranson and Tsimring (2006), was termed motility-induced phase separation. Its theoretical appeal stems from the potential mapping of the nonequilibrium dynamics at large scales to an effective equilibrium theory for the density field Solon et al. (2015a); Nardini et al. (2017); Solon et al. (2018a, b); Prymidis et al. (2016); Paliwal et al. (2017); Bialké et al. (2013); Härtel et al. (2018, 2018). Despite various experimental realizations of self-propelled discs were designed Buttinoni et al. (2013); Deseigne et al. (2010), experimental evidence of this very type of active phase separation is still lacking. Recent experiments with active Janus colloids suggest the hypothesis that polar orientational ordering within clusters interrupts motility-induced phase separation van der Linden et al. (2019).

Orientational symmetry breaking and the emergence of collective motion due to velocity alignment is another central paradigm in active matter – the exploration of phases with orientational order is an integral part of it Vicsek et al. (1995); Toner and Tu (1995); Chaté et al. (2006, 2008); Ginelli et al. (2010); Ramaswamy (2010).11footnotetext: Velocity-alignment may be induced further by hydrodynamic interactions hoell_particle_2018 or dissipative collisions grossman_emergence_2008. We note that torques are typically assumed to be absent for isotropic self-propelled discs. However, torques do only vanish for isotropic, self-propelled particles only if tangential friction among particles and hydrodynamic flows are negligible hoell_particle_2018. The most relevant source of alignment is anisotropic repulsion of spatially extended, self-propelled objects Peruani et al. (2006); Note (1). Recent numerical studies unveiled a large variety of collective phenomena among self-propelled, anisotropic objects including mesoscale-turbulence Wensink and Löwen (2012), formation of bands and aggregates Weitz et al. (2015), accumulation at confining walls Wensink and Löwen (2008) and a complex phase diagram depending on the rigidity Abkenar et al. (2013); Shi and Chaté (2018) and deformability of the particles Menzel and Ohta (2012); Ohta and Yamanaka (2014, 2014); Tarama et al. (2014); Löber et al. (2015). Beyond the inherent theoretical interest in the physics of self-propelled particles, there exists a large number of real-world applications: motile bacteria in two-dimensions Peruani et al. (2012); Wensink et al. (2012); Zhang et al. (2010); Sokolov et al. (2007), biomimetic systems such as motility assays Schaller et al. (2010); Sanchez et al. (2012); Sumino et al. (2012); Tanida et al. (2018); Huber et al. (2018) as well as shaken Narayan et al. (2007); Aranson et al. (2007) and self-propelled, granular rods Kudrolli et al. (2008).

In a nutshell: particle shape controls the physics of active systems – whereas active phase separation is expected for self-propelled discs, anisotropic volume exclusion leads to collectively moving clusters Peruani et al. (2006). Based on that insight, a unified framework encompassing central phenomena observed in active systems can be established as the basis for a theory of active matter. In this context, novel theoretical concepts are called for as the application of concepts from equilibrium statistical mechanics to active matter is limited to a few special cases Solon et al. (2015b); Fodor et al. (2016).

The complexity of models for spatially extended, anisotropic objects has hindered analytical studies and systematic coarse-graining addressing their collective properties – the characterization has been mainly carried out mainly by numerical simulations 222Clustering properties were studied in terms of phenomenological coagulation-fragmentation equations Peruani et al. (2006); Peruani and Bär (2013). The derivation of hydrodynamic equations from microscopic models has only been possible for heuristic point particles with a prescribed velocity-alignment rule based on symmetry considerations and binary collisions Aranson and Tsimring (2005, 2006); Bertin et al. (2006); Ihle (2011); Peshkov et al. (2012, 2014); Großmann et al. (2015); Heidenreich et al. (2016).

Inspired by the quanto-mechanical wave-particle duality, we develop in this work a unifying modeling approach for soft, self-driven agents: each individual entity is represented by an anisotropic, smoothed field whose mutual interactions are derived from the minimization of an overlap energy. Thus, the interaction force and the torque result from a single interaction principle governed by the anisotropy of particles. This approach yields an universal, simple and descriptive model that can be analytically coarse-grained, unifying different paradigms of active matter: motility-induced phase separation and the emergence of orientational order. Combining numerical simulations and analytically derived coarse-grained order parameter equations, we show that motility-induced phase separation cannot emerge for sufficiently anisotropic objects due to the combined action of self-propulsion and anisotropic repulsion. The resulting nonequilibrium stresses acting on the microscale induce orientational alignment of different symmetries locally. In this system, the rigidity of particles determines the symmetry of ordered states: long-lived, giant moving clusters are observed if particles strongly repel each other to prevent particle crossing, whereas large-scale nematic order emerges for particles which can slide over each over. We further indicate that those regimes are separated by a bistable coexistence region, similar to the recently reported ones in motility assays Huber et al. (2018). In short, this study sheds light on the importance of anisotropic repulsion as a source of orientational alignment, particularly on how the interrelation of particle shape, rigidity and self-propulsion determines emergent collective behavior – key elements to be considered in the design of biomimetic materials. Unifying different paradigms of active matter within one theoretical framework is expected to pave the way towards a systematic understanding of soft and deformable active matter such as living cell tissues Löber et al. (2015), bacterial colonies Be'er and Ariel (2019) and driven filaments Huber et al. (2018); Kumar et al. (2018).

II Particle-field representation

Inspired by the wave-particle duality, we represent an active particle at position $\mathbf{r}_{k}$ and orientation $\varphi_{k}$ by an anisotropic, Gaussian distribution $\psi_{k}\!\left(\mathbf{r}\right)$ – see Fig. 1 – as follows:

[TABLE]

This type of representation is reminiscent of the procedure used in smoothed-particle hydrodynamics Gingold and Monaghan (1977); Lucy (1977); Monaghan (2005) as well as Gaussian model potentials to describe molecular interactions Berne and Pechukas (1972). In Eq. (1), the dimensions of the major and minor axis are parametrized by $l_{{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}$ and $l_{\perp}$, respectively, which are oriented along $\mathbf{e}_{\!\!\>{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}\!\left[\varphi\right]=\left(\cos\varphi,\sin\varphi\right)$ and $\mathbf{e}_{\!\!\>\perp}\!\left[\varphi\right]=\left(-\sin\varphi,\cos\varphi\right)$. Throughout, the anisotropy

[TABLE]

or, the aspect ratio $l_{{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}/l_{\perp}$ is the central control parameter for the particle shape determining the collective properties as well. Particles are circular for $\varepsilon\!=\!0$ ($l_{{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}/l_{\perp}=1$) and needle-shaped in the limit $\varepsilon\!\rightarrow\!1$ ($l_{{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}/l_{\perp}\!\rightarrow\!\infty$).

Interaction energy, force and torque

The basic idea is that particles repel each other to minimize their mutual overlap upon encounter Berne and Pechukas (1972); Peruani et al. (2006); Löber et al. (2015). For the Gaussian representation of individual active particles, the overlap $\mathcal{I}_{kj}\!\!=\!\!\int\!d^{2}r\,\psi_{k}\!\left(\mathbf{r}\right)\psi_{j}\!\left(\mathbf{r}\right)$ can be calculated analytically,

[TABLE]

where $\mathcal{Q}=\mathbf{e}_{\!\!\>{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}\otimes\mathbf{e}_{\!\!\>{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}-\mathbf{e}_{\!\!\>\perp}\otimes\mathbf{e}_{\!\!\>\perp}$ (see Supplemental Material (SM) for technical details of the calculation SI ). The interaction energy $\mathcal{U}$ is defined as the sum of binary contributions: $\mathcal{U}=\frac{1}{2}\sum_{k,j}^{N}u_{2}\!\left(\mathbf{r}_{k}-\mathbf{r}_{j};\varphi_{k},\varphi_{j}\right)$, where the binary energy $u_{2}$ is an increasing function of the overlap $\mathcal{I}_{kj}$,

[TABLE]

where $\Delta\mathbf{r}=\mathbf{r}^{\prime}-\mathbf{r}$ is the relative position, $\kappa$ is the interaction strength measured in units of energy and $\mathcal{F}[\xi]$ is a monotonically increasing function of the overlap. This model fulfills two consistency requirements. First, the energy is minimized if particles are nematically (uniaxially) aligned and is therefore in line with Onsager’s mean-field theory Onsager (1949). Second, the energy increases as particles approach each other, hence inducing a repulsive force. In particular, soft and hard objects can be described: if the energy is finite for $\Delta\mathbf{r}\rightarrow 0$, particles are soft whereas these objects can be considered hard if the energy diverges in this limit.

The binary force $\mathbf{f}_{2}\!\left(\mathbf{r},\varphi,\varphi^{\prime}\right)\!=\!-\nabla u_{2}\!\left(\mathbf{r},\varphi,\varphi^{\prime}\right)$ and torque $m_{2}\!\left(\mathbf{r},\varphi,\varphi^{\prime}\right)\!=\!-\partial_{\varphi}u_{2}\!\left(\mathbf{r},\varphi,\varphi^{\prime}\right)$ which a particle at the origin with orientation $\varphi^{\prime}$ exerts on a particle at position $\mathbf{r}$ and orientation $\varphi$ are deduced from the potential energy by differentiation with respect to its position and orientation, respectively (cf. SM SI ). The calculation yields an anisotropic repulsive force between the centers of mass of two particles. Newton’s third law of action-reaction holds for this model. Further, the torque consists of two contributions:

[TABLE]

The first one couples the orientation of a rod to the relative position of another rod: if a rod approaches the center of mass of its interaction partner, it will turn away, thereby minimizing the relative overlap. We refer to this term as collision avoidance. The second contribution is nematic alignment of the body axes Peruani et al. (2008); Ginelli et al. (2010). The prefactors $A$ and $B$ are given in the appendix SI .

We note that forces and torques scale differently in terms of the anisotropy: the force is a zeroth order effect, collision avoidance is a first order effect in the anisotropy $\varepsilon$, whereas nematic alignment scales with $\varepsilon^{2}$ consistently. Accordingly, the torque vanishes identically and the force reduces to an isotropic central body force as expected for spherical particles ($\varepsilon=0$).

Writing force and torque in an analytically closed form (cf. SI ) in contrast to rule-based algorithms is advantageous in several regards. The shape of particles is implicit in the force and torque in this approach. It is thus not necessary to recalculate overlaps or integrate over the cell area in every time step, thereby allowing for a fast numerical implementation. Furthermore, smoothed representations of particles are much easier to handle than hard rods with volume exclusion, reflected by divergent or noncontinuous interaction potentials. Moreover, it enables analytical investigations as discussed in Section IV. Note that the interaction energy is structurally comparable to that used to describe soft deformable particles via a Gaussian interaction potential in Menzel and Ohta (2012); Tarama et al. (2014), however, the resulting expressions differ in details 333To understand in detail the differences, we refer the reader to the respective literature: Eq. (4) in Menzel and Ohta (2012), Eqs. (20)-(24) in Tarama et al. (2014) and Eqs. (7)-(11) in Ohta and Yamanaka (2014) should be compared to the equations given in the Supplemental Material of this work SI ..

Equations of motion

We describe the dynamics of spatially extended, active particles in the overdamped limit as follows Doi and Edwards (1986); Peruani et al. (2006); Baskaran and Marchetti (2008a); Abkenar et al. (2013); Weitz et al. (2015); Shi and Chaté (2018):

[TABLE]

Accordingly, the balance of dissipative and driving force leads to stochastic motion with a mean speed $v_{0}$ in the absence of interactions. The translational ($\hat{\bm{\mu}}$) and rotational ($\mu_{\varphi}$) mobilities are determined by the anisotropy of particles as well as the properties of the surrounding medium; the particular functional forms depend on the experimental context Perrin (1934); *koenig_brownian_1975; Tirado et al. (1984); Doi and Edwards (1986). Similarly, the nature of fluctuations of position and orientation, abbreviated by $\bm{\eta}_{k}(t)$ and $\xi_{k}\!\left(t\right)$, respectively, depends on the context. Here, they are assumed to be Gaussian, unbiased and $\delta$-correlated in time. Furthermore, fluctuations of the center of mass are anisotropic,

[TABLE]

where the diffusion coefficients parallel and perpendicular with respect to the rods’ orientation read $D_{{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}$ and $D_{\perp}$.

We point out that interaction as well as the diffusion matrix and the friction tensor possess “nematic symmetry” – they are invariant under the transformation $\varphi_{j}\rightarrow\varphi_{j}+\pi$ for any rod. Notably, the self-propulsion term $v_{0}\mathbf{e}_{\!\!\>{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}\!\left[\varphi_{j}\right]$ in Eq. (6) breaks this inversion symmetry at the microscale as $\mathbf{e}_{\!\!\>{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}\!\left[\varphi_{j}\right]\rightarrow\mathbf{e}_{\!\!\>{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}\!\left[\varphi_{j}+\pi\right]=-\mathbf{e}_{\!\!\>{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}\!\left[\varphi_{j}\right]$. For this reason, self-propelled rods ($v_{0}>0$) are inherently different from systems without directed self-propulsion on the microscale ($v_{0}=0$), where the diffusive dynamics is invariant under inversions of the orientation vector; this applies to most active nematic models Narayan et al. (2007); Giomi et al. (2011, 2013); DeCamp et al. (2015); Shi and Ma (2013); Thampi et al. (2014); Putzig et al. (2016); Oza and Dunkel (2016); Pismen and Sagués (2017); Cortese et al. (2018).

III Unifying active matter

The smoothed-particle approach allows to unify several seemingly disconnected paradigmatic phenomena reported for different systems before within one framework, such as active phase separation, polarly moving clusters, nematic bands as well as the coexistence of polar and nematic order. Generally, the emergent patterns depend on particle shape and are sensitive to whether particles can slip over each other – implying that the self-propulsion force can overcome the repulsive interaction force – or not Shi and Chaté (2018). The ratio between the strength of the active force and the maximum repulsion force in a binary collision defines the rigidity of the particles. We report the phenomenology of this system in the regime of weak self-propulsion, i.e. high repulsion force implying high rigidity, for different particle shapes in Section III.1 and turn the role played by self-propulsion in Section III.2.

III.1 From self-propelled discs to self-propelled rods

The shape of rigid self-propelled objects determines their collective behavior. Below, we describe the phenomenology observed in numerical simulations along Fig. 2 as the anisotropy of rigid particles is increased.

Breakdown of active phase-separation

In the limiting case of self-propelled discs ($l_{{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}/l_{\perp}=1$), we observe phase separation due to active motion into a disordered, gas-like phase and hexatically ordered aggregates, called motility-induced phase separation (MIPS) Tailleur and Cates (2008); Cates and Tailleur (2015); Nardini et al. (2017); Solon et al. (2018a, b); Prymidis et al. (2016); Paliwal et al. (2017); Bialké et al. (2013); Härtel et al. (2018, 2018). In Fig. 2(a), the enhancement of density fluctuations in the fully phase-separated regime is evident. Since the orientations of discs within an aggregate is disordered, this phenomenon may therefore be described by a scalar field theory for the particle density only Solon et al. (2018b).

Surprisingly, MIPS aggregates are found to melt for small anisotropies of the active particles. In contrast to the phase-separated regime, we observe a drastic decrease in density fluctuations as well as a reduction of local hexatic order already for weakly anisotropic particles (Fig. 2(b), where $l_{{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}/l_{\perp}=1.0424$). See also Appendix A for additional information on the hexatic order.

What is the mechanism behind break up of aggregates? We recall that phase-separation of self-driven spheres arises due to the slow down of particles as they collide. Upon a collision, their orientations point towards the center of clusters such that aggregates are surrounded by a polar boundary layer on average, whereby an active pressure keeps them together [Fig. 3(a)]. Preparing an aggregate in a slap geometry and performing a quench to slightly anisotropic rods [Figs. 3(b,c)] reveals that this polar boundary layer dissolves as a deterministic torque rotates rods away from the boundary of aggregates [Fig. 3(b,c)].

Emergence of orientational order

We observe the formation of states with orientational order by increasing the aspect ratio of particles beyond the breakdown of active phase separation. At first, the system remains globally disordered, however, it displays weak local nematic order as revealed by the semi-analytical analysis discussed in Section IV.0.2. Counterintuitively, local order becomes polar if the anisotropy is increased further, even though the interaction is strictly nematic. The emergence of local polar order is illustrated in Fig. 2(c-e) as well as in Appendix A for intermediate aspect ratios. In particular, polarly ordered, moving domains are found Peruani et al. (2006). Those macroscopic patterns, as shown in Fig. 2(d), are highly dynamic since polar order is inherently related to mass transport, thereby inducing clusters to form, merge and break in a nontrivial fashion Peruani and Bär (2013).

Along with polar domains, topological defects emerge due to collisions of those structures as highlighted by a black dot in Fig. 2(d). Examining the orientation of the rod axis only, i.e. irrespective of the orientation of the self-propulsion force with respect to the body axis, these defects have basically a nematic structure. Defects may, however, be self-motile because of the polarity of directional energy input at the level of individual rods: in Fig. 2(d), a black arrow indicates that rods push towards the center of a $+\frac{1}{2}$-defect thereby creating an active, anisotropic stress, which results in a directed displacement of the defect position; an analogous defect dynamics was observed in Genkin et al. (2017). This mechanism of defect motion in ensembles of polarly driven objects is different from defect motility in active nematics, both dry Narayan et al. (2007) and wet Giomi et al. (2013); Thampi et al. (2014); Pismen and Sagués (2017); Cortese et al. (2018). Furthermore, defects are created and disappear in an intermittent way: due to strong density instabilities, the defects may vanish in the void or penetrate from the boundary of a dense region – the topological charge is therefore not conserved.

Polar domains may become system spanning for intermediate system sizes as shown on the right of Fig. 2 for $l_{{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}/l_{\perp}=4$. These polar bands are comprised of smectic particle arrangements Chaikin and Lubensky (1995); Weitz et al. (2015); Romanczuk et al. (2016). Note, however, that numerical data suggest the absence of long-range orientational order in the thermodynamic limit Weitz et al. (2015).

Let us stress that the observed macroscopic order is polar, while the symmetry of interaction potential is strictly nematic [cf. Eq. (6)]. This implies that the symmetry of the macroscopic order is not imposed by the symmetry of the interaction potential, but emerges spontaneously from the spatial dynamics. Similar behavior was observed experimentally in studies of colonies of myxobacteria as well as corresponding hydrodynamic descriptions Peruani et al. (2012); Harvey et al. (2013).

III.2 Polar vs. nematic order and their coexistence

Let us now fix the aspect ratio of particles and ask what the influence of the rigidity of particles is. In order to answer this question, we varied the self-propulsion force resulting in higher speeds $v_{0}$. In the limit of high self-propulsion, particles may slide over each other upon encounter – an effect reminiscent of tunneling in quantum mechanics which is possible due to the particle-field representation of soft rods – whereas they would be blocked by their interaction partners in the opposite limit. Fig. 4 shows pictorially the phenomenological transition from rigid to soft rods. Large-scale, polar domains are observed for low self-propulsion. In the limit of high activity, in contrast, particles arrange themselves in nematic band-like structures as they are familiar from Vicsek-type, point-like particles with nematic alignment Ginelli et al. (2010). Surprisingly, we found a bistable coexistence region if the order of magnitude of the self-propulsion force is comparable to repulsive forces: nematic bands and polarly ordered domains are observed in an intermittent fashion.

The coexistence of polar and nematic order is detailed in Fig. 4(b). The stochastic switching from polar to nematic states is revealed by anomalously high fluctuations of the polar order parameter. The timescales of these stochastic transitions are remarkable, as they are several orders of magnitude larger than microscopic timescales of the dynamics at the particle level.

The simultaneous existence of polar and nematic states has recently been reported by Huber et al. in Huber et al. (2018) for a motility assay experiment. Those results were rationalized by simulations of self-propelled, flexible filaments which are pulled at one side and interact by volume exclusion. In contrast to this system, the interaction of the self-propelled particles considered in the present study is strictly nematic. The global nematic symmetry is solely broken by the self-propulsion term; individual entities are nematic in contrast to the filaments in Huber et al. (2018). Here, we show for the first time that bistability of polar and nematic structures can also be expected for simple self-propelled rods if the strength of self-propulsion and repulsion are fine-tuned or happens to coincide in a specific application.

We conclude this section by a comment on finite size scaling: we do not expect the bistability of globally polar and nematic states to be retained in the thermodynamic limit as the diffusive motion of particles within the band is too slow to allow for long-range nematic order Weitz et al. (2015). We do rather expect disconnected patches composed of polarly or nematically aligned particles to emerge which may coarsen in a nontrivial way as they interact at the mesoscale. Due to the enormous length-scale separation of several orders of magnitude, a detailed analysis of this phenomenon is beyond the scope of the present work.

IV Analytical insights

The theoretical framework based on a smoothed-particle representation allows for a semi-analytical analysis within kinetic and hydrodynamic theories derived from the microscopic dynamics. In the following, we elucidate the large-scale pattern formation in the course of the transition from active spheres to self-propelled rods in two subsequent steps: the breakdown of MIPS due to combined action of self-propulsion and anisotropic particle shape and the emergence of orientational order from volume exclusion interactions. Afterwards, we conclude by discussing the intricate relation of symmetry of the particle-particle interaction and the symmetries of emergent patterns at large scales in active matter.

The theory is based on the nonlinear Fokker-Planck equation for the one-particle density $P\!\left(\mathbf{r},\varphi,t\right)$,

[TABLE]

with the effective force and torque functionals

[TABLE]

which, in turn, depend on pair distribution function $g_{2}$. We point out at this level that the binary interaction $\mathbf{f}_{2}$ and $m_{2}$, i.e. the force and torque which one particle exerts on another interaction partner [cf. Eq. (5)], are effectively renormalized by the emergent correlations quantified by $g_{2}$. Hence, the interplay of self-propulsion and volume exclusion interactions as well as the collision kinetics can induce novel terms in mesoscale, coarse-grained descriptions, which are not of nematic symmetry, whereas the interaction at the particle level is strictly nematic.

IV.0.1 Breakdown of active phase separation

At first, we study the breakdown of motility-induced phase separation (MIPS) for small particle anisotropies as illustrated in Figs. 2, 3. For the sake of simplicity, we only consider leading orders in $\varepsilon$ (weakly non-spherical limit). This approximation is intended to explain why MIPS aggregates in the classical sense are unstable for anisotropic particles.

In the spherical limit, the force at leading order reduces to an isotropic central body force. Notably, there is a qualitatively novel contribution to the torque at leading order $\varepsilon$ which has not been studied analytically so far in the context of active matter to the best of our knowledge: the torque felt by particle $k$ due to the presence of particle $j$ is proportional to

[TABLE]

cf. Eq. (5). Its meaning and interpretation as collision avoidance is the central subject of the discussion below. Note that the nematic alignment is not relevant in the spherical limit as it is of higher order in $\varepsilon$. Furthermore, mobility and diffusion tensor of individual particles are assumed to be approximately isotropic: $\hat{\bm{\mu}}\approx\bar{\mu}\mathds{1}$ and $D_{0}=D_{{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}\approx D_{\perp}$.

We begin the analysis by calculating the average force and torque – felt by a particle with orientation $\varphi$ –  to first order in gradients:

[TABLE]

The parameters $\zeta_{0,1}$ are nonequilibrium transport coefficients which can be expressed as integrals over the inter-particle correlation function and $\rho$ is the particle density; for details on the derivation, see SI . Numerical measurements of pair-correlations within the particle-based dynamics (see Fig. 8) reveal that the kinetic of particles is such that leads to an enhancement of particle density at the front with respect to the direction of motion of a focal particle: as particles move actively in a semi-dilute environment, they tend to collide with others, and consequently the probability to find a particle in front is higher than in the back. This implies positive transport coefficients $\zeta_{0,1}$. This phenomenon, i.e. accumulation of particle density at the front of the moving particle, was observed in system of self-propelled discs and used to build a scalar field theory to describe phase separation in active disk systems Bialké et al. (2013); Härtel et al. (2018). According to Eq. (10a), the force yields, to first order, an average slowing down of particle speed in high density areas as they bump into their neighbors. The decrease of speed with particle density is the classical mechanism underlying MIPS, first introduced in Tailleur and Cates (2008). In Eq. (10b), we report, however, an important new element not present in classical theories for MIPS: a slight asymmetry of the body shape gives rise to a torque, which induces a rotation away from high density domains.

The physical effects of force and torque become evident at the level of coarse-grained order parameters: the density $\rho\!\left(\mathbf{r},t\right)$ and the polar order parameter field $\mathbf{p}=\left<\mathbf{e}_{\!\!\>{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}\!\left[\varphi\right]\right>$. Equations for the time evolution of these quantities are obtained by performing a mode expansion of the Fokker-Planck equation (8) which yields

[TABLE]

The convective term in Eq. (11a) for the density represents the density-dependent speed reduction due to collisions via $v\!\left(\rho\right)=v_{0}-\bar{\mu}\zeta_{0}\kappa\rho$. A corresponding term is found in Eq. (11b) for the polar order parameter as well. We find that the torque is cast at the field level as an anisotropic, nonlinear flow of the form $\dot{\mathbf{p}}\propto-\varepsilon\Pi_{-}\!\cdot\!\!\;\nabla\rho$ with the tensors $\Pi_{\pm}=\rho\mathds{1}\pm\mathfrak{Q}$, where $\mathfrak{Q}\!\left(\mathbf{r},t\right)$ abbreviates the nematic order parameter field. Accordingly, the torque is such that the density gradients are counteracted by an opposing particle flow. We stress that the polar and nematic order parameter have to be taken into account for the anisotropic, nonequilibrium stresses due to active motion on the hydrodynamic level.

Robustness of MIPS The time-independent solutions of the transport equations (11) imply the polar order parameter $\mathbf{p}$ to be collinear to density gradients $\nabla\rho$, as $\mathbf{p}=D_{0}\!\left(\nabla\rho\right)\!/v\!\left(\rho\right)$ follows from Eq. (11a), cf. the phase-separated state in Fig. 3(a). The theoretical analysis for anisotropic particles reveals, however, that torques will destabilize parallel arrangements of the orientation $\mathbf{e}_{\!\!\>{\mkern 3.0mu\vphantom{\perp}\vrule depth=0.0pt\mkern 2.0mu\vrule depth=0.0pt\mkern 3.0mu}}\!\left[\varphi\right]$ and density gradients $\nabla\rho$ [see Eq. (10b)]. Therefore, this torque, which is proportional to the anisotropy $\varepsilon$, tends to dissolve the polar boundary layer around aggregates. Namely, it induces locally anisotropic stresses whenever density gradients and local order coexist on a coarse-grained level [cf. Eq. (11b)], as we argued before based on numerical simulations (see Section III.1).

To substantiate these arguments, the linear stability of the spatially homogeneous, disordered state is investigated on the basis of Eqs. (11). For self-propelled spheres, the emergence of MIPS is signaled by a long-wavelength instability of the density field Bialké et al. (2013). Here, we examine the stability of the disordered state with respect to long-wavelength fluctuations for anisotropic particles. For this purpose, the dynamics of the polar order parameter field is linearized first by inserting $\rho=\rho_{0}+\delta\rho$ and $\mathbf{p}=\delta\mathbf{p}$. As we are interested in the onset of a long-wavelength instability, the linearized field $\delta\mathbf{p}$ can further be adiabatically eliminated yielding

[TABLE]

In order to determine the onset the long-wavelength instability, it is not necessary to take the fluctuations of the nematic order parameter $\mathfrak{Q}$ into account, which was therefore set to zero. To leading order, one thus obtains an effective diffusion equation $\partial_{t}\delta\rho\simeq\Gamma\Delta\delta\rho$ for the fluctuations of the density around the spatially homogeneous state by inserting this expression into Eq. (11a), where the transport coefficient $\Gamma$ reads

[TABLE]

A long-wavelength instability of the homogeneous state towards a phase-separated state occurs for $\Gamma<0$; the spinodal curves are determined by $\Gamma=0$. Since this instability condition depends explicitly on $\varepsilon$, we can derive a critical particle anisotropy above which MIPS cannot arise (sufficient condition):

[TABLE]

Thus, there is a critical anisotropy beyond which the polar boundary layer which would keep an aggregate together becomes destabilized by torques. In agreement with these theoretical arguments, we find in simulations of the model that MIPS aggregates do not emerge already for rather small anisotropies – for the parameters used in this study, we find the critical values $\epsilon\gtrsim 0.04$, correspondingly $l_{\parallel}/l_{\perp}\gtrsim 1.04$, see Fig. 2(b); see also the dissolution of a MIPS aggregate after a quench from spherical particles to anisotropic rods (Fig. 3). These findings put the relevance of the classical phenomenon of MIPS for self-driven, anisotropic particles, such as self-propelled rods, into question.

IV.0.2 Onset of orientational of order

Let us now examine the emergence of orientational order beyond the spherical limit, as observed numerically for large anisotropies (see Fig. 2). As the interaction at the particle level possesses nematic (front-tail) symmetry [Eq. (4)], one may naively expect the emergence of local nematic order. We will show that the break up of MIPS is indeed followed by a globally disordered phase with local nematic order. Interestingly, local order becomes, counter-intuitively, polar if the aspect ratio is increased even further. In order to identify and understand the emergence of local orientational order at the hydrodynamic level, we derived coarse-grained order parameter equations where hydrodynamic transport coefficients are expressed as integrals over the correlation functions. To simplify the presentation, we concentrate on central, symmetry-breaking terms for the polar and nematic order at the local level, i.e. we expand to lowest order in spatial gradients:

[TABLE]

If the transport coefficient $\sigma_{p}$ is positive, the local polar order parameter grows and, thus, ordered polar structures are expected at local scales. In contrast, the nematic order parameter is relevant at the local level if $\sigma_{p}<0$ and $\sigma_{n}>0$.

The starting point of the analysis is the Fokker-Planck equation [Eq. (8)] for the one-particle density. A mode expansion of this Fokker-Planck equation, where only the relevant terms under consideration are kept, yields the following expressions for the transport coefficients:

[TABLE]

In these equations, we introduced $\tilde{m}_{2}=m_{2}g_{2}$ which is the product of the actual torque $m_{2}$ between two particles and the pair distribution function $g_{2}$. It can be interpreted as an effective mean-field model. We implicitly assumed that the pair-distribution function $g_{2}$ is known and absorbed it into the definition of $\tilde{m}_{2}$. Thereby, the transport coefficients $\sigma_{p,n}$ above still depend on the inter-particle correlations and the collision kinetics.

The simplest mean-field approximation assumes that $g_{2}\approx 1$ – consequently, the effective torque $\tilde{m}_{2}$ is identical to the actual torque $m_{2}$. In this limit, the transport coefficient $\sigma_{p}$ for the polar order parameter is always negative; the integral in Eq. (16a) vanishes for symmetry reasons: $\sigma_{p}=-D_{\varphi}$. The emergence of polar order cannot be described within mean-field theory; it only predicts the existence of an isotropic-nematic transition Baskaran and Marchetti (2008b, a).

How can one rationalize the emergence of polar order provided that the interaction among two anisotropic, elongated particles is strictly nematic? We first recall that only the self-propulsion force breaks the nematic symmetry of the microscopic dynamics [Eq. (6)]. Accordingly, the emergence of polar order has to be related to the collision kinetics and, in particular, the formation of polar clusters Peruani et al. (2006) which is, in turn, reflected by correlations which need to be taken into account properly.

Along with the illustration of polar and anti-polar collisions in Fig. 1, we first give a heuristic argument for the emergence of local polar order in ensembles of rigid rods. Let us consider a collision under an acute angle: an elongated, active particle that collides with a cluster aligns its direction of motion to the local mean orientation. Consequently, it will keep on moving in parallel with this cluster for a significant time. Only rotational diffusion may deflect its direction of motion away from the boundary of the cluster. Therefore, polar clusters are stable and may grow for low rotational noise. In contrast, nematic clusters cannot exist: a particle that collides in an antipolar way with a cluster will just slide off its boundary due to the mobility of the cluster. As the distance of the particle and the cluster increases over time, their positions decorrelate as a result of collisions with other particles. This reasoning leads to the conclusion that the probability to find rigid rods moving in parallel is higher than seeing antipolar arrangements locally as a result of the collision kinetics. This can be verified quantitatively by measuring the pair correlation function numerically, see supplementary Fig. 8(b).

More formally, the effective, binary torque on the field level $\tilde{m}_{2}$, which enters into the relevant hydrodynamic transport coefficient $\sigma_{p}$, is given by the product of the model $m_{2}$ and the correlation function $g_{2}$. In turn, this quantity encodes the kinetics of collisions. Consequently, the effective torque $\tilde{m}_{2}$ may contain new terms which are not present on the particle level: the torque at the particle level contains nematic alignment as $m_{2}\!\propto\sin\left[2\!\left(\varphi^{\prime}-\varphi\right)\right]$, and the pair distribution function contains a positive contribution $g_{2}\propto\cos\!\left(\varphi^{\prime}\!\!\>-\varphi\right)$ as the probability of parallel motion is enhanced. Thus, their product contains effectively positive polar alignment terms proportional to $\sin\!\left(\varphi^{\prime}\!\!\>-\varphi\right)$, contributing to the first Fourier mode $\int_{0}^{2\pi}d\varphi\sin\varphi\,\tilde{m}_{2}\!\left(r,\alpha,\varphi\right)$ in the integral in Eq. (16a), even though polar alignment was not explicitly present at the particle level.

The presence of correlations may renormalize the interaction parameters and even introduce new interaction terms at the field level. Therefore, the coefficient $\sigma_{p}$ can turn positive such that polar terms become relevant on the hydrodynamic level in a model with pure nematic interactions. These arguments crucially depend on the presence of self-propulsion which is the only term that breaks the global nematic symmetry. Accordingly, polar order cannot emerge in the limit $v_{0}\rightarrow 0$.

We investigated numerically how the transport coefficients $\sigma_{p}$ and $\sigma_{n}$ depend on the anisotropy of particles by measuring the pair distribution function $g_{2}$ and evaluating the integrals in Eq. (16). Figure 5 shows the relevant hydrodynamic coefficients together with the variance of the local density as a measure for the fluctuations of the particle density. On the basis of this semi-analytical analysis, we can distinguish the following parameter ranges, which were introduced along with Fig. 2 in Section III.1. Negative values of $\sigma_{p}$ and $\sigma_{n}$ together with a high level of density fluctuations reflects the properties of MIPS as observed for isotropic particles. Increasing the anisotropy, density fluctuations decrease rapidly as MIPS aggregates break up, while $\sigma_{p}$ and $\sigma_{n}$ are negative (local disorder). Subsequently, $\sigma_{n}$ turns positive signaling the emergence of local nematic order, followed by the emergence of local polar order when $\sigma_{p}$ becomes positive.

Thus, a theoretical description in terms of scalar quantities only, such as the particle density, is not applicable for self-propelled rods once vectorial or tensorial order parameters grow at the local level: the polar order parameter field $\mathbf{p}\!\left(\mathbf{r},t\right)$ cannot be adiabatically eliminated if the hydrodynamic transport coefficient $\sigma_{p}$ is positive. In order to capture anisotropic, nonequilibrium stresses on the coarse-grained level – arising due to the combined action of self-propulsion and anisotropic particle shape –, vectorial and tensorial order parameters have to be taken into account. This underlines again the different nature between scalar and anisotropic active matter.

IV.0.3 Discussion: symmetries of emergent ordered states

In Section III, we reported the type of emergent patterns to be essentially dependent on the spatial dynamics, namely whether rods can possibly slip over each other or not. To rationalize this observation, we consider the spatial dynamics [Eq. (6a)] comparing the order of magnitude of the active force propelling a particle forward with the characteristic magnitude of the repulsive force which is determined by the ratio of an energy scale $\kappa$ and the typical size of a rod. If the spatial dynamics is dominated by the active force, rods may easily slide over each other on a timescale which is faster than the timescale of the rotational dynamics, implying that the orientational dynamics is slow. Therefore, one may simplify Eq. (6a) to $\dot{\mathbf{r}}_{k}\simeq v_{0}\mathbf{e}\!\left[\varphi_{k}\right]$. In this limit, it turns out that the large-scale behaviour of spatially extended rods becomes similar to the dynamics of Vicsek-like rods with nematic alignment Ginelli et al. (2010); Shi and Chaté (2018). Hence, the following phenomenology is expected: for low rotational angular diffusion, a spatially homogeneous, nematic phase emerges at the mesoscale (finite system size). By increasing the noise, the level of nematic order decreases. Close to the order-disorder transition, the system demixes into a high density region which is nematically ordered and a low-noise area where particles move in a disordered fashion. This is the type of pattern shown on the right of Fig. 4. We therefore conclude that the nematic alignment term in the torque [Eq. (5)] dominates the large-scale dynamics in the limit of high activity. Accordingly, corresponding field theories for Vicsek-type self-propelled rods Peshkov et al. (2012); Großmann et al. (2016) may account for the observed pattern formation phenomena, such as band formation and the emergence of nematic order. Apparently, the fact that particles push each other is of minor importance in this parameter regime. In other words, positional correlations are less relevant if particles move fast as the system becomes well-mixed, i.e. particle positions decorrelate quickly when rods can slip over each other.

By decreasing the self-propulsion force or, equivalently, increasing the repulsion strength, particles would get blocked upon encounter, positional correlations build up and the pair-correlation function becomes increasingly relevant such that mean-field arguments are not applicable. In this regime, patterns may emerge whose symmetry differs from the symmetries of the microscopic interaction. It remains an open challenge for future work beyond the present study to derive the pair correlation function of anisotropic, self-propelled objects from first principles – analogues to corresponding theories for self-propelled discs Härtel et al. (2018) – in order to explain the emergent bistability in parameter regimes where self-propulsion and repulsion are comparable which is most difficult to assess analytically as standard series expansions fail.

V Summary & outlook

We introduced a new approach to study active systems: inspired by the quanto-mechanical wave-particle duality, active particles are represented by smooth fields and particle-particle interactions are derived from the minimization of an overlap energy. Given the nature of the model both force and torque can be analytically obtained. Importantly, this modeling technique allows studying the transition from self-propelled discs to active elongated, anisotropic, soft objects by performing continuous deformations on the shape as well as on the rigidity of the active particles. Here, we analyzed analytically how motility-induced phase separated phases for isotropic, active particles become unstable for weakly anisotropic, self-propelled objects. Specifically, we demonstrated that the combined action of anisotropic repulsion and self-propulsion leads to the emergence of effective torque terms which – above a critical particle aspect ratio – dissolve the polar boundary layer required to maintain motility-induced phase separated aggregates. These findings provide an understanding of the role played by particle anisotropy regarding the robustness of the active phase separation described in terms of scalar field theories for the particle density (MIPS).

Furthermore, we used the same theoretical framework to study the onset of (orientationally) ordered phases. We established that both particle aspect ratio and rigidity – in comparison with the self-propulsion – control the symmetry of the pair correlation function. We demonstrated that the emerging order, either nematic or polar, is not dictated exclusively by the symmetry of the interacting potential but rather given by the resulting symmetry of both, interaction potential and the pair correlation function. That is why polar and nematic structures can simultaneously coexist in a system of identical particles with a purely nematic interaction potential. This can be realized by fixing the rigidity and varying the particle aspect ratio, or fixing the aspect ratio and varying the rigidity. At the level of the analytically derived coarse-grained order parameter equations, this would require showing that it is possible to observe nematic order and vanishing polar order as well as non-vanishing polar and nematic order simultaneously. In the latter scenario, the presence of nonlinear cross-coupling terms in the evolution of nematic and polar order parameters Peshkov et al. (2012) implies that the magnitudes of these coefficients control whether nematic or polar orders prevail, or both coexist. All these situations are possible for the studied system of anisotropic, self-propelled particles. In short, our analysis reveals that the symmetry of macroscopic order is an emergent and dynamic property of active systems, as recently suggested from the analysis of a motility assay experiment Huber et al. (2018).

The developed framework allows studying motility-induced phase separation for isotropic, self-propelled particles, its breakup with particle anisotropy as well as the counterintuitive emergence of both, polar and nematic macroscopic order for the same type of active particles, thereby providing a unified picture of most relevant phenomena reported in active matter systems. We therefore expect our framework to shed light on a large number of applications, including collective dynamics in collection of cells, the growth of bacterial colonies and self-organized patterns in systems of deformable, active filaments. Natural extensions of the developed approach range from the addition of hydrodynamic flows to the study of polydisperse systems, among many others.

Acknowledgements.

RG and FP acknowledge support from the Agence Nationale de la Recherche via Grant No. ANR-15-CE30-0002-01. RG acknowledges funding from the People Programme (Marie Curie Actions) of the European Union’s Seventh Framework Programme (FP7/2007-2013) under REA grant agreement n. PCOFUND-GA-2013-609102, through the PRESTIGE programme coordinated by Campus France. ISA was supported by the NSF 446 PHY-1707900.

Appendix A Additional Langevin simulations

Snapshots

In addition to Fig. 2, the gradual emergence of local polar order for intermediate values of the anisotropy is illustrated in Fig. 6. In both cases, the transport coefficient $\sigma_{p}$, which indicates the growth of local polar [cf. Eq. (15)] order, is positive.

The breakdown of a phase-separated state (MIPS) for anisotropic particles as shown in Fig. 2 is accompanied with the breakdown of local hexatic order. Figure 7 shows the same states like Fig. 2, but the color coding corresponds to the local hexatic order parameter.

Pair distribution functions

A central theoretical argument is based on the enhancement of the probability to find a particle in front (with respect to the direction of self-propulsion) as compared to the probability to find it in its back. This is reflected by the pair distribution function $g_{2}$, specifically by the Fourier component

[TABLE]

as shown in Fig. 8 (a) for spherical particles, compare Bialké et al. (2013); Härtel et al. (2018).

In Section IV.0.2, we argued further that the probability of parallel motion is enhanced on average as a consequence of the anisotropic body shape Peruani et al. (2006), see also the schematic collisions in Fig. 1. This is also reflected by the pair correlation function. In particular, a positive contribution to the Fourier component

[TABLE]

reflects that the probability of moving together in groups ($\varphi\approx\varphi^{\prime}$) is larger than swimming in an anti-parallel fashion ($\varphi\approx\varphi+\pi$). This argument has been verified by numerical measurements of the respective part of the pair distribution function, cf. Fig. 8(b).

Movies of Langevin simulations

We provide three movies SI for the nontrivial states shown in Fig. 2: mov1.mp4 (2a), mov2.mp4 (2d) and mov3.mp4 (2e). Furthermore, the animation mov4.mp4 corresponds to the simulation which is depicted in Fig. 4(b).