Tuning Self Organization Of Confined Active Particles By Steepness Of The Trap

TL;DR

This study explores how the steepness of a confining trap influences the collective behavior and ordering of self-propelling particles, revealing a transition from boundary flocking to spatially ordered structures as the trap becomes less steep.

Contribution

It demonstrates how tuning the trap steepness induces a transition from boundary flocking to hexagonal spatial order in confined active particles.

Findings

Particles form boundary flocking in steep traps.

Decreasing trap steepness leads to hexagonal order.

Order stability persists after trap removal.

Abstract

We consider collective dynamics of self-propelling particles in two dimensions. They can align themselves according to the direction of propulsion of their neighbours, together with a random perturbation (i.e. rotational fluctuation). They are also interacting with each other by a soft, isotropic, excluded-volume interaction. Particles are confined in a circular trap. The steepness of the trap is tuneable. Their average packing fraction and strength of rotational fluctuation are low. When the trap is steep, particles flock along its boundary. They form polar cluster that spreads over the boundary. The cluster has no spatial or structural order. We show, when the steepness is decreased beyond a threshold value, the clusters become round, compact and eventually spatial order (hexagonal) emerges in addition to the pre-established polar order within them. We investigate kinetics of such ordering. We find that while rotating around the centre of the trap along its circular boundary, the clusters need to roll around their centre of mass, to be spatially ordered. We also discuss the stability of the order when the trap is suddenly switched off.

Equations

$$U(r_i)=U_0(\tanh (q(r_0-r_i))+1)$$

$${\hat{\mathbf{n}}}_i=R_{\xi }(\theta )\circ \left(\frac{1}{N_{S_i}}\sum _{k\in S_i}{\hat{\mathbf{n}}}_k\right)$$

$$V_{wca}$$

$$m\frac{d{\mathbf{v}}_i}{dt}=-\gamma {\mathbf{v}}_i-{\mathbf{F}}_i+f_0{\hat{\mathbf{n}}}_i+\sqrt{2\gamma K_{B}T}{\eta }_i$$

$$\text{Pe}=\frac{f_0}{\gamma \sqrt{K_{B}T}}.$$

$$\Lambda =\frac{f_0}{\gamma \sqrt{ϵ}}$$

$$S=q\sigma .$$

$${\psi }_{6,i}=\frac{1}{N_i}\sum _j\exp (\mathbf{i}6{\theta }_{ij})$$

$$\mathcal{P}=\left|\frac{1}{N}\sum _i{\hat{\mathbf{n}}}_i\right|$$

Full text

Active colloids, Confinement, Phase transition

††preprint: APS/123-QED

I Introduction

Active (self-propelled) interacting particles are prone to be phase separated to form clusters with varied dynamical and structural properties, even in two dimensions (2D). By now this has been extensively studied in bulk via theory and experiments Deseigne et al. (2010); Cates et al. (2010); Peruani (2012); Schwarz-Linek et al. (2012); Palacci et al. (2013); Liu et al. (2017); Theurkauff et al. (2012); Buttinoni et al. (2013); Redner et al. (2013); Cates and Tailleur (2013); Fily and Marchetti (2012); Farrell et al. (2012); Peruani et al. (2011, 2010); Ginot et al. (2018); Peruani et al. (2006); Tung et al. (2016); Weitz et al. (2015); Stenhammar et al. (2014); Fily et al. (2014); Levis and Berthier (2014); Stenhammar et al. (2013); Mognetti et al. (2013); Patch et al. (2017); Martín-Gómez et al. (2018); van der Linden et al. (2019); Sese-Sansa et al. (2018); Shi and Chaté (2018); Alarcón et al. (2017). Passive equilibrated systems lack self-propulsion and self-aligning interactions, particularly of non-steric origin. Hence the dynamics and the patterns observed in the active systems are usually more intricate than the systems in equilibrium. Though surprisingly, it turns out that intriguing physics of two step melting (solid to hexatic to liquid) Halperin and Nelson (1978); Nelson and Halperin (1979); Strandburg (1988); Bernard and Krauth (2011); Engel et al. (2013); Kapfer and Krauth (2015), a characteristic feature of 2D passive (equilibrated) systems, may be important to understand the dynamics of clustering and the structures within, observed in self-propelling but not self-aligning soft discs Klamser et al. (2018a, b). Recent simulations suggest that, similar to equilibrium systems, self-propelling discs, while interacting with each other by soft, isotropic, repulsive potentials can exhibit hexatic phase between solid and liquid phases while melting, even for considerably high activity. According to Klamser et al. (2018a), the system needs to be highly active to open up another criticality (independent of the melting one) where self-propelling particles, thanks to their motility, come together to phase separate and form clusters coexisting with gas-like fluid phase. The motility induced phase separation (MIPS) Cates and Tailleur (2015) is the predominant feature of the active particles where self-propulsion and isotropic excluded volume interactions dominate over self-aligning interactions. It is an outcome of a positive feedback between particle-accumulation and their motilityStenhammar et al. (2013); Gonnella et al. (2015); Cates and Tailleur (2015); Wittkowski et al. (2014); Marenduzzo (2016). Though MIPS is reported frequently in the current literature but its characterisation remains illusive. Various issues related to MIPS, e.g. whether the dense phase obtained in MIPS is liquid-like Klamser et al. (2018a) or solid-like Bialké et al. (2012); Redner et al. (2013) or polycrystalline Cugliandolo et al. (2017), whether the kinetics of phase separation in all the active systems are equivalent or not Redner et al. (2013); Buttinoni et al. (2013); Briand and Dauchot (2016) etc., are debatable.

The other extreme, where aligning interaction dominates over excluded volume interaction in a system of self-propelling particles, with sufficiently low fluctuations and high density, transitions from isotropic to ordered phases are observed, instead of MIPS. Depending on the types of the aligning interactions Vicsek et al. (1995); Chaté et al. (2006); Liebchen and Levis (2017), choice of the parameter space Caussin et al. (2014), underlying symmetry of the system together with non-equilibrium fluctuations at large scales Solon and Tailleur (2013); Solon et al. (2015a); Solon and Tailleur (2015) various kinds of ordered phases (e.g. - nematic, polar, smectic) and ordered structures (e.g. - periodic density waves, solitonic bands, polar-liquid droplet, rotating micro-flocks etc.) Solon et al. (2015b); Liebchen and Levis (2017) are obtained. The transition from isotropic to polar phase among a 2D populations of self-aligning agents with no excluded volume interactions, is observed in the seminal works, e.g. Vicsek et al. (1995); Toner and Tu (1995). The more complete understanding regarding the transition came much later Grégoire and Chaté (2004). Recently we know that the transition from spatially homogeneous isotropic phase to spatially homogeneous polar phase occurs via a phase coexistence with microphase separation (which gives rise to an effective smectic order) Solon et al. (2015a). There are several studies where the ordered phases and their properties, such as giant number fluctuations are explored in experiments Nishiguchi et al. (2017); Narayan et al. (2007).

Often in the experiments it is hard to disentangle the effect of steric and aligning interactions. This motivates one to analyze the collective dynamics and spatiotemporal pattern formation within active systems in presence of both. Various complex dynamical patterns, e.g. polar clusters of varied sizes, band-like or lane-like motile, amorphous structures, flowing crystalline clusters etc. can emerge as a result of the competition between self-propulsion, excluded volume effect and self-alignment interactions Martín-Gómez et al. (2018); Shi and Chaté (2018). It has been claimed in Sese-Sansa et al. (2018) that velocity alignment in such active systems can promote MIPS by reducing the effect of rotational diffusion. Though according to the studies in Weber et al. (2014); Shi and Chaté (2018), structural order can also be reduced to a smaller region of the full phase diagram due to aligning interactions. In recent experiments, in one hand we have evidence showing that aligning interaction can interrupt full phase separation via MIPS van der Linden et al. (2019), on the other hand it has been shown that a polar flock with excluded volume interaction can undergo a freezing transition in bulk being assisted by MIPS Delphine et al. (2019).

There are remarkable effects, observed within the collective dynamics of self-propelling units, where the confining boundary plays a major role. For example, accumulation of active particles close to the boundary, often leads to non-Boltzmann density profiles which in some simpler cases can be calculated exactly Tailleur and Cates (2009, 2008). It signifies the out-of-equilibrium character of the system. When an active particle faces a boundary, it cannot turn back from it unless its velocity points away from the boundary. Therefore they are often stuck to it. If there are many such particles, they can come to the boundary one after the other and accumulate. This makes the turning back from the boundary even harder for the particles in front. This effect is observed in varied set-ups, for example — Bricard et al. (2015); Vladescu et al. (2014); Wensink and Löwen (2008); Elgeti and Gompper (2009); Volpe et al. (2011); Wysocki et al. (2015). The steric interactions between the boundary of the trap and the particles are enough for their accumulation at the boundary without the need for hydrodynamic interaction Elgeti and Gompper (2013). Though the hydrodynamic flows close to the boundary are essential to establish crystalline order within the accumulation resulting from flow induced phase separation (FIPS) Singh and Adhikari (2016); Thutupalli et al. (2018). The stability of the order depends on the geometry of the flow fields which essentially depends on the boundary conditions. Inter-particle as well as particle-wall hydrodynamic interaction both can be important for phase separation in a dense, active, confined colloidal systems Zöttl and Stark (2014). Apart from crystallisation, in the presence of confinements, other complex patterns, such as active fluid pumps and vortices are also formed within various confined, simulated active systems (e.g. active Brownian particles and run and tumble particles Nash et al. (2010); Hennes et al. (2014); Menzel et al. (2016)) and in real active systems (e.g. swimming bacteria, colloidal rollers, vibrated polar discs etc.Wioland et al. (2013); Lushi et al. (2014); Bricard et al. (2015); Kudrolli et al. (2008); Deseigne et al. (2012).

Here we focus on a subtle but simple way to introduce structural order (hexagonal) in addition to the polar order, already established within a motile cluster of confined, self-propelling, soft, circular discs interacting with each other by velocity alignment interaction. No long ranged interaction of hydrodynamic origin is assumed among the particles. It may be noted here that in many cases hydrodynamic interactions do not produce any qualitative difference in the self organisation of active particles (e.g.-Kaiser et al. (2017) ).The structural order is established without affecting the motility and polarity of the cluster. The mean packing fraction of the system is kept as low as 15% of the total surface area. The trap to confine the particles and the inter-particle steric repulsion both are spherically symmetric. In order to develop the structural order within the cluster we will exploit the interplay among their self-propulsion, aligning interactions and boundary effects, in particular, the steepness of the trap boundary.

II Main results

Before going into the details of the model and the results, in this section we will state our main findings briefly. We explore the influence of the confinement on the collective dynamics of active particles to control their organisation. In general, the interplay among activity, self-aligning interaction and steric repulsion offers complex, collective dynamics with pattern forming, non-equilibrium steady states in bulk. Here our focus will be on the role played by the confining boundaries. We demonstrate the influence of the confinement with a model representing soft, active discs with self-aligning interaction and confined within a circular trap. Structural (hexagonal) as well as directional (polar), both order can be established in the system. The phase diagram of the system with respect to both, in the parameter space of activity, boundary properties for a fixed strength of the fluctuations in self-aligning interaction are studied in detail. Through out the study we have kept the packing fraction low ($\simeq$15 % of the surface area) and constant. We find that the fluctuations having strength beyond a threshold value, can destroy both the order. More importantly, we also find that below the threshold, though the polar order is stable with respect to the steepness of the trap, but the structural order is unstable with respect to that.

We demonstrate this by simulating the system from same random initial configuration but with two different steepness of the circular trap (everything else e.g., depth of the trapping potential, packing fraction, activity, self-aligning and steric interactions etc., are kept same for the two), which finally produces two differently organised steady states — for low steepness the final state is structurally ordered and for high steepness there is no structural order (Fig[1]). In both the cases particles accumulate close to the boundary and form dynamic, polar clusters that rotate along the boundary of the circular trap. The size of the clusters of the active particles depends on the length scales set by the trap. To emphasise the influence of the steepness of the trap boundary further and to achieve a control on the degree of the structural order, one can also start from the polar but structurally disordered phase and reduce the steepness slowly (quasistatically) over time to develop the order in addition to the polar order within the cluster. This can be used as a strategy to control the degree of the order within the polar cluster of self-propelling and self-aligning active particles in a confinement. In Fig[2] we have shown that how the structural order varies with the steepness of the trap.

After identifying the relevant parameter spaces for various phases in the phase diagrams we explore the kinetics of the ordering. We find that the transition from the isotropic to the polar phase in our system resembles the same within a population of Vicsek agents. The transition is robust even in the presence of the confinement and steric interactions. With decreasing steepness of the trap, we find that the structural order can emerge within the self-organised polar phase of the active particles, keeping its polar order intact. The structural order is the outcome of the interplay among the spatially non-uniform repulsion from the boundary, self-propulsion and self-aligning interactions present in the system. When the boundary of the trap is very steep, particles within the clusters are rotating along the boundary. The interface between the cluster and the trap is smooth. At any point of time, the velocity direction of individual interfacial particle points almost along the tangent of the circular trap at that point. With reducing steepness, while rotating, interfacial particles start to bounce at the boundary of the trap repeatedly. Consequently the smooth, circular interface turns into rough by creating ripples. With further decrease of the steepness, the ripples grow and the polar cluster is fragmented into smaller ones. Smooth fronts or interfaces being unstable, can break into clusters due to different types of dynamical instabilities. For example, in a suspension of colloidal micro-rollers fingering instability, originated from purely long-ranged hydrodynamic interactions, can break smooth shock front into small clusters Driscoll et al. (2017). In the system of our concern, after fragmentation, the smaller clusters also bounce repeatedly at the trap boundary with low steepness. Consequently a force couple acts within the clusters to make individual cluster rolling while rotating along the boundary. Due to rolling, particles within a cluster are facing the confining force from all directions. Eventually the clusters become compact enough to support the structural order, keeping its polarity and collective motion intact. These polar and structurally ordered clusters can merge together to become a single cluster or can collide with each other. Though these events can drastically destroy the orders within a cluster, but it is a transient phenomena. After one such event, if a cluster is allowed to relax to a steady state without facing further collision or merging events, eventually it regains both structural and polar order. It implies that the orders, which are developed within the clusters, are robust against such events. Hence we obtain a simple protocol to control the degree of structural order within a polar cluster which can be implemented in experiments. Pertinent issues regarding crystalline order within confined active particles are raised and investigated in Briand et al. (2018). There, initially self-propelling and spontaneously aligning hard discs are allowed to form a close-packed configuration with perfect hexagonal order within a hexagonal arena. Removing particles from the configurations causes condensation of shear along the stacking faults, that leads to large scale shear flow within the system. Surprisingly both structural and polar order can survive against the flow, providing us flowing, polar crystal. But the packing fraction used in the study is very high ( $\simeq$81%) compared to the study here. Consequently, we will see that the kinetics of ordering here is fundamentally different from Briand et al. (2018). We have also explored the stability of the structural order when the confinement is suddenly turned off and showed that the decay of the order can follow power law over time. The power can vary over time intervals, indicating multiple regimes of the decay, crossing over one by the other. The details of the results will be provided later.

III Model

We consider a collection of $N$ self-propelling particles at $\{{\bf{r}}_{i}\}$ with velocities $\{{\bf{v}}_{i}\}$ where $i=1,2,...,N-1,N$. They are confined within a two dimensional circular trap modelled as a continuous, spherically symmetric function of $r_{i}=|{\bf{r}}_{i}|$ as,

[TABLE]

such that $2U_{0}\geq U(r_{i})\geq 0$. We consider that the boundary of the trap is at $r_{0}$, where the potential felt by a particle is : $U(r_{0})=U_{0}$. The energy of the individual self-propelling particle remains always lower than $U_{0}$ such that they can never cross the circular barrier. The parameter $q$ controls the steepness of the boundary. The extent of the repulsion from the boundary towards the centre of the trap can be quantified by the length scale $l=q^{-1}$. It is small if $q$ is large and vice versa. Though the magnitude of the repulsive force from the boundary will increase as $q$ becomes larger. Therefore, large $q$ limit is close to the hard wall limit where the layer of the particles closest to the boundary faces a very high repulsive force but the layers afterwards face a very small repulsion directly from the wall. On the other hand, when $q$ is small, the potential gradient and consequently the force gradient due to the boundary become smoother, making it soft. In case of the soft boundary, the extent of the repulsive force from the boundary is larger. In this case, not only the layer of the particles closest to the boundary but the layers behind can also feel the repulsion from it. Though its magnitude decreases gradually to zero as one goes towards the centre of the trap from its boundary. In this work, we focus on how the gradient which is controlled by a single parameter $q$ of the potential, affects the self-organisation of self-propelling particles confined within it. Through out the work we consider $qr_{0}>>1$ such that the influence of the boundary on the particles is considerable only when they are close to it and become negligible as we go far from it i.e. in the central region of the trap.

While propelling with constant force $f_{0}$, an individual particle try to follow its neighbours within a certain circular area of radius $r_{int}<<r_{0}$ surrounding the particle. In order to do that the particle follows Vicsek-like aligning interactionVicsek et al. (1995); Weber et al. (2014). The particle decides its direction of active force ${\bf{\hat{n}}}_{i}$ to be aligned with the average velocity direction of its neighbours including itself. This local alignment is not perfect due to the presence of an angular noise $\theta$ randomly chosen between [math] and ${2\xi}{\pi}$ at every time step. Here $\xi$ is a dimensionless parameter that lies between [math] and $1$ and defines the strength of the angular noise. One can express the alignment interaction as,

[TABLE]

where $\hat{\bf{n}}_{k}={\bf{v}}_{k}/{|{\bf{v}}_{k}|}$, $S_{i}$ is the neighbouring circular area of $i$th self-propelling particle, within which there are $N_{S_{i}}$ particles that include the neighbours of $i$th particle and $i$th particle itself. They are indexed with $k=1,2,....,N_{S_{i}}$. $R_{\xi}(\theta)$ is the operator that rotates (denoted by $\circ$) the average direction of the neighbours within $S_{i}$ by the randomly chosen angle $\theta$ to model the imperfection in the alignment. Instead above, one may also consider continuum version of Vicsek dynamics Farrell et al. (2012); Liebchen and Levis (2017); Martín-Gómez et al. (2018) though the qualitative features of the system will remain same Chepizhko and Kulinskii (2010).

Within a distance smaller than $r_{int}$, the particles are also repelling each other by a spherically symmetric force derived from a soft, continuous potential $V_{wca}$: Weeks-Chandler-Andersen (WCA) interaction Weeks (1971) to account for the effect of excluded volume interaction between the particles,

[TABLE]

Here $r_{ij}$ is the distance between i-th and j-th particles and $\epsilon$ determines the strength of the repulsive potential between them. The length scale involved in the repulsive potential is $\sigma$ which is the measure of the size of an individual particle. We have considered $\sigma$ as the unit of length. We consider $r_{int}=2\sigma$ to ensure alignment even when particles are very close to each other. Here we assume that the alignment is of non-steric origin which may be realised by the systems having spherically symmetric vibrated disks Scholz and Pöschel (2017) or active colloids under magnetic or long-ranged hydrodynamic interactions respectively. It can also be mentioned here that the model described here captures the essential features of active systems such as migratory cells Stichel et al. (2017), where both the effects of excluded volume and alignment interactions are relevant.

At time $t$, if the position and velocity of $i$th particle are denoted by ${\bf{r}}_{i}(t)$ and ${\bf{v}}_{i}(t)$, its dynamics can be written within Langevin paradigm as,

[TABLE]

where $m$ is the mass of a particle, $\gamma$ is the friction coefficient due to the surrounding fluid and ${\bf{\eta}}_{i}$ is the thermal noise, i.e. $\langle\eta_{i}\rangle=0,\langle\eta_{i}(t_{1})\eta_{j}(t_{2})\rangle=\delta_{ij}\delta(t_{1}-t_{2})$. We consider mass of the particle to be unity and $\gamma$ to be fixed at high value ($\sim 10$) such that for unit $\epsilon$ and $\sigma$, $\sqrt{\epsilon}/\gamma\sigma<<1$ which is motivated by the systems with high viscosity (e.g. colloidal suspensions, bacterial suspensions etc.). The force ${\bf{F}}_{i}$ is the superposition of all the forces on $i$th particle, provided by the trapping potential $U$ and collisions with other particles via WCA interaction, i.e. ${\bf{F}}_{i}=-\nabla_{i}U+\sum_{j}{\bf{F}}_{ij}^{wca}$. $f_{0}$ is the magnitude of self-propelling force, the direction of which is decided by Eq.[2]. The strength of the thermal noise is related to the friction coefficient and temperature of the surrounding fluid by fluctuation-dissipation relation which ensures the equilibrium of the system in absence of self-propulsion. We consider the system to be far from equilibrium due to self propulsion where the effect of passive, thermal diffusion is small. Therefore through out the calculation we consider $K_{B}T/f_{0}\sigma\leq 10^{-4}$.

In similar systems, role of the defects on the verge of the transition from fluid like disordered phase to ordered phases (for example, hexagonal and polar phases) have been explored in bulk Weber et al. (2014). Recently, while varying different parameters, it has been shown that similar model can also accommodate other complex steady states in bulk, e.g. finite-sized polar clusters, band-like or lane-like highly dynamic structures and their hybrids Martín-Gómez et al. (2018). Collective dynamics of self-propelling chiral entities are investigated by similar models with Vicsek-like alignment combined with active torque Liebchen and Levis (2017). Here we will explore the collective dynamics of such self-propelling, mutually aligning, finite-sized particles under the influence of a circular confinement. We will show that the steepness of the trap has major impact on the self organisation of such active particles. In particular, for a fixed packing fraction of confined active particles, depending on the Peclet number and the strength of the rotational noise, by tuning steepness of the trap one can induce structural disorder-order transition into them. We will keep the packing fraction $\phi={N\sigma^{2}}/{r_{0}^{2}}=0.15$ fixed through out the work to investigate the effects of alignment, excluded volume and confinement in detail.

At this stage it is convenient to identify few non-dimensional parameters that control the collective behaviour of the confined particles. Individual particle is experiencing forces from other particles due to collisions which is of the order of ${\epsilon}/{\sigma}$; it is experiencing thermal force of the order of ${K_{B}T}/{\sigma}$ due to the surrounding fluid at temperature $T$; then from the boundary of the trap it is facing repulsive inward force of the order of $qU_{0}$. Finally due to active degrees of freedom, individual particle exerts force $f_{0}$ to their surrounding. To ensure confinement, here we consider that ${qU_{0}}/{f_{0}},{qU_{0}\sigma}/{\epsilon}$, ${qU_{0}\sigma}/{K_{B}T})\geq 1$. Again, we consider ${\epsilon}/{K_{B}T},{f_{0}\sigma}/{K_{B}T}$ to be larger than unity such that we can focus on the interplay between steric interaction and activity only.

For the system here, we define Peclet number as,

[TABLE]

which quantifies the self-propulsion with respect to thermal diffusion by taking the ratio of friction limited active velocity and thermal velocity of the particles. With $\text{Pe}\rightarrow 0$ the system becomes passive (equilibrated with temperature $T$), losing their active degrees of freedom ($\{{\bf{\hat{n}}}^{t}_{i}\}$) and related dynamics.

In the system of our interest, there are competing aligning and de-aligning interactions. The competition can be quantified by their respective ranges. The ranges of Vicsek-like aligning and simultaneously de-aligning interaction are given by $s_{1}=\pi r_{int}^{2}$ and $s_{2}=\xi\pi r_{int}^{2}$ respectively. The competition between them can be quantified by the parameter $s_{2}/s_{1}=\xi$ that varies between [math] and $1$. At $\xi\rightarrow 1$ limit the particles can be maximally de-aligned due to rotational fluctuations in the alignment interaction and for $\xi\rightarrow 0$ the opposite happens. One can note here that for alignment interaction to work, particles should see each other, i.e. the average distance between them should be at least $r_{int}$ or less. As we keep the packing fraction low ($\phi=0.15$), they can see each other only after a certain critical value of Pe. We will be back to this point while discussing the phase diagrams in later section. The limit where $\xi\rightarrow 1$, the self-propelling particles resemble active Brownian particles (ABP) Fily and Marchetti (2012); Cates and Tailleur (2015); Redner et al. (2013); Stenhammar et al. (2013); Stenhammar (2014); Buttinoni et al. (2013) confined within the trap. In 2D, ABP encounter MIPS in bulk at much higher density as compared to here. In the regime of density and Pe we are concerned here, as ABP lack alignment interactions, it will not produce clusters beyond the critical size required for MIPS.

The particles here can behave as effectively soft as well as hard, depending on certain condition. The inter-particle steric interaction between any two particles here becomes increasingly repulsive as they come close and it diverges when their centres overlap i.e. at $r_{ij}=0$. Therefore, the particles can be soft only up to a certain extent and beyond that it becomes extremely hard such that the particles cannot come closer to each other. The hardness of the particles can be parametrised as

[TABLE]

If $\Lambda>1$ the particles are soft and if $\Lambda\leq 1$, it is hard. Here, varying $f_{0}$ while keeping the other parameters fixed, we have explored hard as well as soft particles. For very large $\Lambda$ and Pe (in the limit of $\epsilon\rightarrow 0,K_{B}T\rightarrow 0$), the particles are close to Vicsek agents. In bulk, the homogeneous state of Vicsek agents can be transformed to a polar state when $\xi$ is lower than a critical strength $\xi^{c}$ via a first-order phase transition and when trapped, they can exhibit different complex dynamical patterns (e.g. swarms) Vásárhelyi et al. (2018).

Finally, the dimensionless steepness parameter of the circular wall of the cavity can be written as

[TABLE]

We will show that one can tune $S$ by varying $q$ to control the self-organisation of the particles externally.

In the next section onwards, after mentioning simulation details, we will present the phase behaviour of the system, with respect to activity (Pe) and the steepness parameter ($S$) of the trap for fixed value of $\xi$ that supports polar order. Next we will discuss the kinetics of phase transformations together with the stability of the ordered and polar phases and finally we will conclude by summarising the work together with its possible experimental realisation.

We simulate the system by integrating the equations of motion given in Eq[2] and Eq[4] simultaneously with open boundary condition. To update velocity and position we follow velocity Verlet algorithm with $dt=10^{-3}$. All the data presented here are tested with total number of particles $N=5\times 10^{3}$ and $N=10\times 10^{3}$ particles to estimate the finite size effect which has turned out not to cause any qualitative change to the data presented here. To analyse the data we have used OVITO Stukowski (2009).

IV Results

IV.1 Order Parameters

We begin by analysing the steady states of the system obtained by varying Pe, and $S$ at fixed packing fraction of the particles. The steady states of the system are characterised by their spatial and directional organisations. The spatial organisation of the particles is quantified by the local (at the level of individual particle) hexagonal order parameter, defined as

[TABLE]

where ${\bf{i}}=\sqrt{-1}$ and $N_{i}$ is the number of topologically (Voronoi) nearest neighbours of $i$th particle. The summation is over all such neighbours of $i$th particle and $\theta_{ij}$ is the bond angle between $i$th and $j$th particles with respect to some arbitrarily chosen reference frame. For a perfect triangular lattice in 2D, $Q_{6}^{i}=|\psi_{6,i}|$ is unity. With disorder its value reduces from unity and for a maximally disordered configuration it is zero. For a given configuration of particles, with $Q_{6}^{i}$ one can quantify local order and disorder by identifying ordered zones and defects in between. To quantify the global positional order of the system in a steady state we consider global hexatic order parameter as $Q_{6}=\langle|\langle\psi_{6,i}\rangle_{i}|\rangle_{\text{(r,t)}}$ where $\langle.\rangle_{i}$ denotes the averaging over all particles and $\langle.\rangle_{(r,t)}$ denotes averaging over realisations and time in steady state. The state of the system having $Q_{6}$ above (below) a certain threshold is ordered (disordered).

The local polar order of individual particle can be quantified by the angle made by $\hat{\bf{n}}_{i}$ with respect to an arbitrarily chosen reference frame and to quantify the global polar order, following Vicsek et al. (1995); Weber et al. (2014); Martín-Gómez et al. (2018) we define the polar order parameter as,

[TABLE]

where ${\mathcal{P}}\equiv|{\bf{P}}|=|P_{r}{\hat{r}}+P_{\theta}{\hat{\theta}}|=\sqrt{P_{r}^{2}+P_{\theta}^{2}}$ in plane polar coordinate. When all the particles are propelling being perfectly aligned with each other, $\mathcal{P}=1$ and when they are moving randomly, $\mathcal{P}=0$. In a steady state, if the time and realisation averaged value of $\mathcal{P}$ goes beyond a certain threshold, the state is polar and otherwise non-polar.

The polar state here is the reminiscent of the flocking state of Vicsek agents Vicsek et al. (1995), which can be stable when the value of $\xi$ is below a threshold value. In this state, when the particles encounter the boundary of the trap, the radially outward component of the force exerted by the particles to the trap is balanced by radially inward repulsion from the boundary. Also in the polar state, as $\xi$ is low, it is more probable for a particle to remain close to the boundary than to turn around from it Tailleur and Cates (2009). Thus the particles in the polar state accumulate close to the boundary and flock together along it, giving rise to self-rotating clusters in the trap. The sense of rotation is spontaneously chosen. In a similar set up where hydrodynamic and electrostatic interactions are also taken into account, formation of other complex structures such as vortices are also reported Bricard et al. (2015). They are not finite size effect but genuine, macroscopic, non-equilibrium phases of active matter in confinement Caussin (2015). The major concern of our work is the emergence of structural (e.g. hexagonal) and polar order within such dynamic clusters.

Due to circular symmetry of the boundary, in a highly polar steady state, we find that $P_{r}<<P_{\theta}\simeq{\mathcal{P}}\rightarrow 1$. In the state where the particles are moving randomly, both $P_{r}$ and $P_{\theta}$ are very small, close to zero, causing ${\mathcal{P}}\rightarrow 0$. Therefore in the system of our concern, emergence of polarity can equivalently be quantified by $P_{\theta}$ alone as $P_{r}$ is negligibly small in both polar and non-polar states.

IV.2 Phase Behaviour

We analyse the phase behaviour of the system by calculating $Q_{6}$ and $\mathcal{P}$ in steady states with different values of Pe and $S$. We fix $\xi=1/2\pi$. In steady states we calculate $Q_{6}$ to estimate the structural order that varies along $\xi=1/2\pi$ plane of a 3D phase diagram having three axis: Pe, $S\phantom{x}\text{and}\phantom{x}\xi$. This variation is plotted in Fig.[3A]. Similarly in steady states we also calculate $\mathcal{P}$ to estimate the polar order that varies along $\xi=1/2\pi$ plane Fig.[3 B]. Below we will discuss them one after the other.

(a) According to Fig.[3A], we will discuss various phases with varying $S$ but fixing different ranges of Pe. For all $S$ but at low Pe (range : 1 $\leq$ Pe $\leq$ 8 ), we find that the particles are randomly distributed through out the trap. In this regime, thermal fluctuations dominate over activity, prohibiting accumulation close to the boundary. The value of $Q_{6}$ is quite low ($Q_{6}\leq 0.4$).

As we raise Pe (range : 8 $\leq$Pe$\leq$ 15), self-propulsion starts to dominate over thermal fluctuations and the particles accumulate close to the boundary forming clusters. The corresponding value of $Q_{6}$ becomes higher than the earlier. It remains between 0.4 and 0.55. The degree of ordering is still low. There is a remarkable difference between the shape of the clusters at $S$ lower than a threshold value $S_{c}$ ($S_{c}\simeq 0.70$) and clusters at $S>S_{c}$. At higher steepness the clusters wet the circular wall of the trap i.e. the contact angle between the layer of the particles and the boundary of the trap is very small De Gennes et al. (2013) . On the other hand, if the steepness is reduced, the clusters become somewhat round, increasing the contact angle between the clusters and the wall of the trap.

Next, when Pe is increased further (range : 15 $\leq$ Pe $\leq$ 40), we find a clear order-disorder transition, together with the shape transformation of the clusters as mentioned before, when we vary $S$ beyond a threshold value $S_{c}$. In this range of Pe, the steady states for $S<S_{c}$ ( $S_{c}\simeq 1.0$), the clusters are round in shape and they have polycrystalline order with dynamic grain boundaries Godrèche et al. (1991). The average value of $Q_{6}$ in this regime (around $S\rightarrow S_{c}^{-}$) can become as high as 0.85. While increasing $S$, beyond $S_{c}$, the clusters become strip-like, wetting the boundary of the trap. They loose their structural order. Consequently, $Q_{6}$ decreases up to 0.55.

Next, Pe is increased even further (range : 40 $\leq$Pe$\leq$ 200). In this regime there is always an order-disorder transition if one goes beyond $S=S_{c}$ while increasing $S$. Interestingly we find that in this range, as Pe becomes gradually higher, $S_{c}$ becomes gradually lower. At Pe$=$200, the transition happens around $S_{c}\simeq 0.5$ whereas at Pe$=40$ it was taking place around $S_{c}=1.0$. From lowest value of S considered here is $0.1$. In this range of Pe, from the lowest value of $S$ up to $S=S_{c}$ we obtain round-shaped cluster with hexagonal order as earlier and beyond we obtain strip like clusters along the circular wall with reduced structural order ($Q_{6}\simeq 0.55$). This indicates that for very low and very high values of Pe, degree of the structural order within the clusters is low. Only the values of Pe in between (range : 20 $\leq$Pe$\leq$ 200) can facilitate the ordering if the steepness of the trap is below $S_{c}$.

(b) From the phase behaviour with $\mathcal{P}$ (Fig.3B) we find that two distinct regions separated by a boundary. For Pe$>10$ and $S>0.2$ we find the system is polar ($\mathcal{P}\simeq 0.9$). To reach to the polar phase, the particles have to $\it{see}$ each other. As the packing fraction is low, the Pe of the particles required to interact with each other should be beyond a threshold which is around $10$ here. The polar clusters rotate along the circular boundary of the trap with spontaneously chosen direction (clockwise/anticlockwise). On the other hand for Pe$<10$ the particles are moving randomly and therefore we obtain isotropic phase with $\mathcal{P}\leq 0.2$. In this regime thermal fluctuation plays a dominant role which prevents clustering. One can note here that the transition from isotropic to polar phase is almost independent of the steepness of the trap, unless $S<0.2$. For very low steepness the transition happens at higher Pe.

IV.3 Kinetics of Ordering

From now onwards we will focus on how the structural order emerges in the system with Pe $\geq 20$ and when $S$ is decreased from the values larger than $S_{c}$ up to the values less than $S_{c}$, keeping $\xi$ fixed at $1/2\pi$. The value of the packing fraction $\phi$ is also fixed at $0.15$. At this parameter space the system becomes polar. From the phase diagram before it is apparent that even if the steepness of the trap is varied over a broad range, it is not able destabilise the polar order rather while lowering, it introduces another order—the structural (hexagonal) order—in addition to the polar order. The onset of such order is analysed below.

It is apparent from Fig3(A) that there is a threshold value for the steepness ($S_{c}$) of the trap below which the structural order can emerge for Pe $\geq 20$. When $S>S_{c}$, apart from very low Pe, for all other Pe, particles accumulate at the boundary of the trap, spreads over it and flock together along the boundary. Thus we obtain polar clusters rotating along the circular boundary of the trap. The boundary of the trap, when steep (i.e. $S>S_{c}$), can exert repulsive force only on the interfacial particles i.e. the peripheral particles of the cluster which are closest to the trap boundary. The particles behind the interfacial particles are flocking together along the circular boundary under active force, inter particle collisions and thermal fluctuations. The direct influence of the repulsion from the boundary on them is very small. It is apparent from the Fig4(A,a), where we have shown the position of the particles of such clusters with their individual direction of self propulsion together with the values of the potential they are facing from the confining trap. When $S\geq S_{c}$, the particles follow the direction almost tangential to the trap. We find the interface between the cluster and the trap smooth Fig4(A,a).

When the steepness of the trap approaches $S_{c}$, ripples start to appear at the interface, making it rough Fig4(B,b). With lowering $S$, together with the interfacial particles, the other particles behind the interface start to feel the radially inward repulsion from the boundary. The repulsion decreases exponentially over a length scale $l$ as one goes radially inward, from the boundary of the trap. Therefore the particles feel spatially non-uniform repulsion from the boundary: particles closer to (far from) the boundary face more (less) repulsion from it. This causes the particles which are close to the boundary at a certain instant of time, to go away from it in the next instants, where the repulsion is very small. But from there, due to self propulsion they again come close to the boundary and be repelled. Therefore due to the interplay between the self propulsion and non-uniform repulsion from the boundary of the trap, the particles are bouncing back and forth at the boundary. As the cluster is polar, even while bouncing, the particles are following each other efficiently. Therefore they bounce coherently at the boundary of the trap. Consequently ripples are formed at the interface of the polar cluster and the circular boundary of the trap. When $S$ is decreased further, more and more particles face non-uniform, radially inward repulsion from the trap boundary. Size of the ripples grows. Roughness of the interface enhances that finally leads to fragmentation of a cluster into smaller ones. This phenomena is clearly demonstrated in Fig4. The fragmented clusters are also bouncing back and forth at the boundary due to non-uniform repulsion from the boundary and self propulsion of the particles. It is apparent from the oscillation of the position of their centre of mass over time, shown in Fig5(A, inset). In the same figure we have also plotted the position of the centre of mass of the clusters that we have when the boundary of the trap is very steep i.e. $S>S_{c}$. The oscillation is negligibly small comparison to the former which implies the absence of bounce in the latter.

When $S<S_{c}$, we have found that the fragmented polar clusters are repeatedly bouncing at the boundary of the trap. Therefore the radial component of the velocities ($v_{r}^{i}$) of the particles within the clusters can be positive as well as negative. The distribution of $v_{r}^{i}$ over the particles of the clusters are shown in Fig5(A, B, C). It shows that the particles of a cluster approaching the wall of the trap have positive radial velocities, particles at the wall have zero radial velocities and particles bouncing off the wall have negative radial velocities. Most of the particles in such cluster have non-zero $v_{r}^{i}$. Almost half of the particles of the cluster has $v_{r}^{i}\geq 0$ and the other half has $v_{r}^{i}<0$. Therefore the distribution $P(v_{r}^{i})$ (shown in Fig5(A)) has the mean at zero and its standard deviation $\sigma_{v_{r}}$ is very large. On the other hand, most of the particles of the clusters in the trap with steep boundary ($S>S_{c}$) have zero or very low radial velocities. The distribution $P(v_{r}^{i})$ in this case, has zero mean and very small $\sigma_{v_{r}}$ in comparison to the former.

For the particles in a highly viscous suspensions, the velocity is proportional to the force where the proportionality constant is the mobility $\mu=\gamma^{-1}$ of the particles. Therefore the distribution of the radial component of the net force (say, $\mathcal{F}_{r}^{i}$) on the particles should be qualitatively same as radial velocity distribution. This implies that for almost half of the cluster that formed within a smooth ($S<S_{c}$) trap, $\mathcal{F}_{r}^{i}\geq 0$. Particles of this half are approaching the boundary before the bounce. For the other half $\mathcal{F}_{r}^{i}<0$. These particles are going away from the boundary after the bounce. Consequently a force couple is generated within the bouncing cluster Fig5(D). It generates a torque to roll around its own centre of mass. The force couple is absent or negligibly small within the clusters confined by the steep trap where $S>S_{c}$ because most of the particles have zero or very small radial velocities. Therefore no torque is generated within such clusters to roll.

One can find the evidence of the rolling of an individual cluster while plotting a typical trajectory of any particle within the cluster (Fig5(F)). The trajectory on XY plane manifests superposition of two time periodic functions with different frequencies and phases. One of them represents its rotation (flocking) along the boundary which is about the origin of the trap. The other one represents its rolling about its own centre of mass. It is expected that this feature is absent for the cluster within a trap having $S>S_{c}$ (Fig5(E)) as $\sigma_{v_{r}}$ is small in that case. In Fig5(G), for steep as well as smooth traps, we have calculated the frequency of rolling ($\omega_{i}$) with respect to centre of mass for the individual particles and show how they are distributed. From there it can be inferred that for steep traps, the variation of $\omega_{i}$ over the particles in the cluster (can be quantified by the standard deviation $\sigma_{\omega}$) is much smaller than its variation observed in the clusters within smooth traps. It implies that for most of the particles in a cluster within a steep trap, the rolling frequency is negligibly small. Only the interfacial particles have finite rolling frequency which is not enough for the whole cluster to roll. On the other hand, in case of the clusters within a smooth trap, almost all the particles have finite rolling frequency ($\sigma_{\omega}$ is large). Thus the whole cluster rolls around its centre of mass while rotating around the centre of the trap. Rolling makes the cluster round in shape.

One can also find how $\sigma_{v_{r}}$, for which the force couple and thereby the rolling of a cluster occurs, scales with the steepness of the trap. The energy involved in rolling of a cluster, $e_{r}\propto I\omega_{c}^{2}$ where $I$ is its moment of inertia and $\omega_{c}$ is its rolling frequency. The width of a cluster measured from its interface with the trap boundary depends on the steepness of the trap. It is proportional to $l=q^{-1}$ which implies, $I\propto m_{c}l^{2}$ where $m_{c}$ is the mass of the cluster. So $e_{r}\propto m_{c}l^{2}\omega_{c}^{2}$. Again, from previous discussion we know that the rolling of a cluster occurs due to a force couple. At moderately high Pe, when the clusters are formed, the magnitude of the active force $f_{0}$ exerted by individual particle, contributes to the force couple much more than other forces present in the system. In the scale of a single cluster, $f_{0}$ can be considered as the magnitude of active force per unit area. Therefore the total active force involved in a cluster is proportional to $f_{0}l^{2}$. Thus the energy involved in rolling can also be written as, $e_{r}\propto f_{0}l^{3}$. Equating the expressions for energies involved in rolling we obtain $\omega_{c}\propto f_{0}^{1/2}l^{1/2}=f_{0}^{1/2}q^{-1/2}\propto f_{0}^{1/2}S^{-1/2}$. From previous discussion, we know that the cluster rolls because of the large variation of radial velocities of the particles i.e. for large $\sigma_{v_{r}}$. Therefore, for simplicity if one assume $\sigma_{v_{r}}\propto\omega_{c}$, it provides $\sigma_{v_{r}}\propto S^{-1/2}$ for constant $f_{0}$. We have tested the scaling relation between $\sigma_{v_{r}}$ and the steepness of the trap $S$ and the result is in Fig5(H).

From the discussions above, it is apparent that for $S<S_{c}$ and moderately high Pe, the rotating (along the circular boundary) polar clusters can roll around their centre of mass. We have shown that how the rolling frequency of a cluster scales with the steepness $S$ of the trap. The rolling is important. We find in our simulations that if a cluster rolls, then only structural order emerges within it. This is because the structural order to emerge with short range repulsive forces such as WCA, the particles should be close enough with each other making the cluster compact. In the current context, the circular boundary can compactify the rotating clusters only when they can roll. When $S>S_{c}$, though the clusters rotate along the boundary of the trap they cannot roll around its centre of mass. Therefore the particles within the cluster, which are away from the boundary of the trap, never face the confining potential directly. So they cannot be close enough to each other to exhibit the order. On the other hand, when $S<S_{c}$, for moderately high Pe, the clusters can rotate as well as roll along the boundary, facing the confinement through out its periphery over time. Thus it can be compactified by the circular boundary of the trap to manifest the order.

V Stability of order

In the sections before, we have discussed that tuning the steepness of the trap can give rise to the hexagonal order within the population of the self-propelling particles, in addition to the pre-established polar order. The relevant parameter space where the transition can occur, is identified and the kinetics of the onset of such transition is also detailed. Below we will discuss stability of the hexagonal order in absence of the trap.

Self-propelling but non-aligning particles can be trapped acoustically and when the concentration is high they crystallise to a close pack structure within the trap. It has been demonstrated that in such cases when the trap is switched off, the crystal explodes and homogeneous distribution of active particles are obtained quickly, within 12 seconds after switching of the trap Takatori et al. (2016). Here the hexagonal order is the outcome of the interplay between the self-propulsion and velocity alignment interactions of the particles together with the steepness of the trap. The ordering takes place at low average packing fraction. Such an order is unstable if any one of the trio is absent. One can investigate that how the order being unstable decays over time, when the trap is suddenly switched off. We consider that initially we prepare the system such that the particles are in a hexagonally ordered, polar cluster. Then suddenly the trap is switched off and we calculate $Q_{6}$ over time. In Fig6 we plot the time evolution of $Q_{6}$ after switching off the trap. We find that $Q_{6}$ decays with power law. Interestingly we find that there are two distinct regimes of the decay where it follows two different powers. The regimes crosses over each other. It indicates that possibly multiple physically different processes are responsible for the decay. More systematic study is needed along this line to identify the exponents and the physical processes, which is not included in the present article.

VI Conclusion

The model and the results discussed here is amenable to the experiments. For example, one may consider induced-charge, electrophoretic, self-propelled Janus colloidal particles, as described in van der Linden et al. (2019) within a confinement providing electrostatic repulsion to the particles. The repulsion decays exponentially with the distance between the particle and the confining boundary. The decay involves a length scale (Debye length) that typically depends on the solution, the size of the particles and the zeta potential of the surface of the colloidal particle and the boundary Rashidi and Wirth (2017). Thus the decay length may be tuned to control the steepness of the confining boundary, which is essential for experimental realisation of the system discussed here.

We have considered a circular (thereby, spherically symmetric) trap with tuneable steepness, which can confine self-propelling, self-aligning, soft discs at low packing fraction. The particles being active, accumulate at the boundary of the steep trap. They form strip-like polar cluster that spreads and rotate along the boundary. The polar order is stable because their rotational fluctuation is quite small compared to self propulsion and velocity aligning interaction. Tuning the steepness of the trap we have shown that as the steepness goes below than a threshold value, the cluster becomes round in shape and hexagonal order emerges within it, keeping the polar order intact. The transition occurs because, as the steepness of the trap is reduced, the polar, strip-like cluster is fragmented into smaller clusters which not only rotate but also roll along the boundary of the trap. It is the rolling for which the trap can push the particles inside a cluster, from every direction, as it rotates along the trap boundary. Due to this, in addition to the polar order, the cluster becomes compact enough to exhibit hexagonal order in the presence of short-ranged, isotropic, inter-particle repulsion. The hexagonal order is unstable in the absence of the boundary. When the trap is switched off, the order decays over time with power law.

VII ACKNOWLEDGMENTS

AS thanks University Grants Commission Faculty Recharge Program (UGCFRP), India.