Spontaneous flows and dynamics of full-integer topological defects in polar active matter

Polar active matter of self-propelled particles sustain spontaneous flows through the full-integer topological defects. We study theoretically the incompressible flow profiles around (cid:2) 1 defects induced by polar and dipolar active forces. We show that dipolar forces induce vortical flows around the +1 defect, while the flow around the (cid:3) 1 defect has an 8-fold rotational symmetry. The vortical flow changes its chirality near the +1 defect core in the absence of the friction with a substrate. We show analytically that the flow induced by polar active forces is vortical near the +1 defect and is 4-fold symmetric near the (cid:3) 1 defect, while it becomes uniform in the far-field. For a pair of oppositely charged defects, this polar flow contributes to a mutual interaction force that depends only on the orientation of the defect pair relative to the background polarization, and that enhances defect pair annihilation. This is in contradiction with the eﬀect of dipolar active forces which decay inversely proportional with the defect separation distance. As such, our analyses reveals a long-ranged mechanism for the pairwise interaction between topological defects in polar active matter.


Introduction
Active matter refers to non-equilibrium systems of interacting, self-propelled entities that consume energy from their surrounding in the form of persistent motion and their collective interactions lead to emergent, dynamical patterns, and self-sustained flows 1,2 .Models of active matter are largely inspired by biological systems from bacterial suspensions 1,2 and cell monolayers 3,4 and down to subcellular active systems such as mixtures of cytoskeletal filaments and motor proteins [5][6][7][8] .However, this also pertains to non-living active systems such as layers of vibrated granular matter, microrobots, or synthetic catalytic nanomotors 2,9 .
Several hydrodynamic models have been proposed to capture See DOI: 10.1039/cXsm00000x/ the macroscopic dynamics and emergent phenomena corresponding to a collection of active (self-propelled) particles with different symmetries and alignment interactions 2,[10][11][12] .The prototypical models are based on the analogy to liquid crystals formed by rod-like particles with polar or nematic symmetries in their alignment interactions.Active rods with only orientational alignment act as headless "shakers" and form nematic phases described by a slowly-varying director field n which has head-tail (nematic) symmetry 1,2,13,14 .By contrast, active rods that align their direction of motion tend to flock into polar systems that are described instead by a slowly-varying polar vector field p [10][11][12]15 . Actve rod-like particles generate persistent flows sustained by the active stress originating from extensile/contractile dipolar forces.On hydrodynamic scales, this active stress is proportional to the nematic Q tensor order parameter, i.e. σ a = α 0 Q with a proportionality coefficient α 0 as an effective activity parameter 2,16 .Topological defects, innate to ordered systems with broken continuous symmetries, are also present in active systems.The interplay between active stresses and nematic distortions feeds into self-sustain flows and the proliferation of topological defects to generate chaotic flows also known as active turbulence 17? ,18 .
From the rotational symmetry modulo π of the Q tensor, it follows that the lowest energy orientational defects have a ±1/2 topological charge corresponding to ±π jump in the orientational phase around them.The +1/2 defects acquire a self-propulsion due to the net active flow passing through their cores 19,20 .The stable, low-energy defects in polar active systems have instead ±1 topological charges corresponding to a 2π phase jumps around the defects.This is analogous to vortices in the XY-model of 2D ferromagnets 15,21 .
The bulk of recent studies have focused on the formation and characterization of half-integer nematic defects (see recent reviews 22,23 ).This is in part due to the ubiquitous emergence of the nematic defects in a wide range of biological systems from subcellular filaments 8,24 to bacterial colonies 25,26 and assemblies of eukaryotic cells 27,28 .This is despite the fact that several biological active entities, such as bacteria or eukaryotic cells, are endowed with a clear head-tail asymmetry and directional selfpropulsion, which characterizes a polar order for these systems.Such a polar order at the scale of collective is apparent from flocking domains within bacterial colonies 26 and eukaryotic cells 29 .Nevertheless, because of the appearance of half-integer defects, at the coarse-grained level, these systems are often modelled as active nematics neglecting the polarity of the self-propelled particles.There have been proposed models that couple the evolution of polar and nematic order parameters [30][31][32] to allow for the coexistence of both types of symmetries.Recently, in Ref. 33 , a hydrodynamic model was proposed for the coexistence of both nematic and polar alignment interactions in the polarization field p, and was used to study different active turbulence regimes sustained by both half-and full-integer topological defects.In addition to dipolar (nematic) active forces, the model also includes a polar active force f = αp p, which describes the self-propulsion direction and is shown to suppress defect-laden active turbulence 34 .
More recently, theoretical and experimental studies have revealed the importance of full-integer topological defects in cell assemblies.In particular, it is shown that positive full-integer defects formed due to collision of two nematic half-integer defects in fast-moving bacterial colonies can lead to the verticalization of bacteria and escape to the out-of-plane directions 26,35 .Furthermore, positive full-integer defects induced by confinement of myoblast cells in circular geometries are shown to activate cell differentiation and formation of 3D helical structures 36 .A corresponding theoretical analyses, in the limit of compressible flows inside the core region of defects in small confinements, have revealed the corresponding flow fields and force patterns, and shows how confinement-induced topological defects can be used to probe the material properties of the cell layers 37 .
In this paper, we carry out a theoretical analysis that reveals subtle cross-talks between polar and dipolar active forces in generating spontaneous incompressible flows around both ±1 defects.We theoretically predict spontaneous vortical flows induced by dipolar active forces around +1 defects which correspond to isotropic active stress and pressure fields.Interestingly, the competition between dipolar force and viscous force leads to a flow reversal close to the core of the +1 defect.This effect is also confirmed by numerical simulations of the full hydrodynamic model including additional passive stresses.However, it turns out the hydrodynamic screening induced by friction lifts up this flow reversal effect.By contrast, the −1 defects which have an innate 4-fold symmetry in the polarity field lead to 8-fold symmetries of the active flows induced by dipolar active forces.The same 8fold symmetry is present also in the profile of the pressure field.We demonstrate that polar active force trigger a distinct active flow pattern characterised by uniform flow in the far-field of its source, i.e. ±1 defects.This is very different that the active flows sustained by dipolar forces which decay algebraically with distance.Polar active forces have a drastic effect on the mutual interactions between defects by promoting fast annihilation of defect pairs through non-local and non-reciprocal attraction forces.We show that this distance-independent mechanism of defect pair annihilation is responsible for the suppression of active turbulence by polar forces as recently reported in Ref. 34 .
The paper is organized as follows: We begin in Section 2 by introducing the flow equations within a minimal polar hydrodynamic model and derive the corresponding defect kinematic equations using Halperin-Mazenko formalism 38 .The main analytical results on the active flow velocity induced around isolated ±1 are discussed in Section 3. We also compare the analytical predictions with direct numerical simulations and find very good agreement of the flow profiles around defects and the flow reversal pattern near the +1 defect cores.In Section 4, we discuss the polar active force and its effect on the defect kinematics.In particular, we consider the motion of a pair of oppositely charged full-integer defects under the polar active flows induced by each defect and demonstrate that polar active forces enhance the defect annihilation rate.In a recent study Ref. 34 , it was numerically evidenced that polar active forces suppress defect-laden active turbulence and tend to restore polar order.Here, we demonstrate theoretically that polar active forces have a net effect on the defect kinematics to promote defect pair binding and subsequent annihilation of defects of opposite topological charges.The final section provides a summary of the theoretical insights and concluding remarks.

Hydrodynamic model of polar active matter
We consider the hydrodynamic model 33 that describes the collective dynamics of self-propelled entities in terms of the evolution of the polar order parameter p given by where γ is the rotational viscosity, λ is the flow-aligning parameter 39 , and the free energy favoring the polar order is described as Here, Kp is the isotropic elastic constant for distortions in the polarity field and A is the height of local energy barrier.The polarity is coupled with the flow field u which is described by the incompressible Stokes equations where the dipolar active forces are mediated by the nematic order through the nematic tensor Q ij = (p i p j − p 2 2 δ ij ) and the activity parameter α 0 .The incompressibility constraint determines the fluid pressure P .The viscosity is set by η while the friction with a substrate is introduced by the frictional drag Γ 40 .Since we focus on theoretical derivations of the active flows induced by polar and dipolar forces, we hereby neglect the additional passive stresses which depend on the polarity and that are typically present in the Stokes equations simulated numerically.
We consider the dimensionless forms of these equations following the appropriate rescalings in units of length ξ = Kp/A and time τ = γ/A.The system is then controlled by three dimensionless parameters: the scale number ζ = d /ξ is the ratio between the hydrodynamic dissipation length d = η/Γ and the nematic coherence length ξ, and the rescaled activity α0 = α 0 τ /(ξ 2 ) = α 0 γ/(Kp).Thus, the dimensionless flow equations read as where the dipolar active force induced by the active stress is For infinite systems, the flow velocity and pressure originated from dipolar active forces can be calculated from convolution integrals of the source terms and the corresponding Green's functions as and To derive analytical expressions, we evaluate the source term due to the active stress using the order parameter for an isolated point defect located at the origin and given by where θ = arctan(y/x) is the polar angle in a coordinate system centered at the defect position and χ(r) is the core function, where we assume that the core size is much smaller than any other length scales and set χ(r) = 1.The angle of the polar vector order parameter for an ideal defect is given by Θ(θ) = qθ + φ 21 , where qθ is the singular part giving the winding number q when performing the integral dΘ = 2πq on a contour surrounding the defect and φ is a constant.For defects of charge q = +1 this constant sets the orientation of the defect and can be ignored since we can always transform to a system with φ = 0 by a change of basis.When q = +1 it is impossible to remove a non-zero φ by changing the basis 41 .The positive defects are rotationally symmetric and the constant φ becomes important since it distinguishes different types of defects: φ = 0 corresponds to an aster and φ = π/2 gives rise to a vortex, and any value in-between corresponds to a spiral defect.For a more intuitive depiction of this, we plot in Fig. 1 a) and c) the p field on the circumference of a circle centered at ±1 defects to show that φ corresponds to the constant angle that p makes with the radial direction.For the negative defect one can define a polarisation from the line where p is pointing radially outward from the core, the angle of this polarisation is half of φ.The different patterns of the p field around a ±1 defect are also illustrated in Fig. 1 b) and d).
We use the Halperin-Mazenko formalism 42,43 for tracking topological defects as zeros of the polar order parameter to derive the corresponding equations of motion of defects from the evolution of the p field similar the approach from Ref. 38 to describe orientational defects in active nematic films.The basic idea is that since the ±1 defects are associated simultaneously with topological singularities in the orientation field and zero magnitude of the vector order parameter, we can track their position by Dirac delta functions centered at the zeros of the p field.Hence, a configuration of well-separated defects punctuating the p-field corresponds to a defect charge density field, which can be written equivalently either as a superposition of the delta functions associated with the topological singularities located in the physical space or as the zeros in the order parameter space, where D is the determinant of the polarity distortion tensor ∇ p, i.e.D = ij ∂ i px∂ j py which can be expressed equivalently in the complex representation ψ = px + ipy as 2iD = ij ∂ i ψ∂ j ψ.The D-field is zero in regions of uniform polar order and becomes nonzero where there are distortions in the orientation field.For configurations of well-separated defects punctuating a uniform polar order state, the D field is zero everywhere except at the defect positions labeled by the index n, where the topological charge qn = ±1 is determined by the sign of D. Thus, the D field represents a non-singular charge density field.
It can be shown that the D-field follows the conservation law 42 with the current density J = ij (∂ t ψ∂ j ψ) determined by the evolution of the polar order.Combining this with the conservation of topological charge density ρ, we can derive a general expression for the defect velocity in terms of the D and its current, vn = J/D| rn 38 .We parameterise the polarity field from the viewpoint of the n-th defect as ψ = ψne iφn , where φn is the smooth phase at the defect position rn and ψn = | r − rn|e iθn using that the polarization vanishes linearly while its phase θn is singular at the pointwise defect position 44,45 .Within this approach, the general expression of the defect current density can be reduced to a closed expression for the defect velocity given by Using the stationary phase approximation, i.e. that the constant (equilibrium) phase remains stationary when it is punctu- ated by moving defect singularities, we can further simplify the defect velocity and express it in terms of the net spontaneous flow and the forces induced by the other defects as 38,46 vn = u( rn) The topological defects interact through Coulomb-like forces, where like-signed defects repel, and opposite-signed defects attract each other.However, there are additional interactions through the flow field u which depends on the dipolar active forces.In the next section, we derive analytic expressions for this flow velocity and discuss its effect on the defect motion.

Active flow fields around ±1 topological defects
From the parameterization of the p-field for a pointwise defect, we can evaluate the dipolar active force induced by an isolated defect with charge q = ±1 in an otherwise homogeneous polarity field with constant background orientation φ as where r ⊥ = (y, −x), and R 2φ is a matrix that rotates the vector by 2φ, and can be removed by a change of basis.Notice that the first term in F + a is a source of gradient flow, which however is removed by pressure through the incompressibility condition.The second term proportional to sin(2φ) is related to a rotated gradient and induces a purely vortical flow.
We insert these forces into the integrals in Eq. ( 8) and ( 9), and solve them as described in the Supplementary Material 47 .The resulting expressions for the active flow velocity and pressure for a +1 defect reduce to Here L is a cutoff scale set by the system size.We have written the velocity field as a complex field defined as u = ux + iuy, where r and θ are the polar coordinates centered at the defects position.We define the scaled radial coordinate as r = r/ζ, which is the same as changing the length scale to the hydrodynamic dissipation length d (see Fig. 2, with the cutoff scale L = 50).Notice that the pressure is isotropic and its gradient force cancels the radial component of the dipolar active force, thus no net pressure flow.Furthermore, the corresponding vorticity is also isotropic and is non-zero at the centre of the defect as evidenced also in Fig. 2.This has important consequences for the stability of vortex, spiral, or aster shaped positive defects as discussed in detail in Ref. 15 .The sign of the global circulation is modulated by the character of the +1 defect through φ.
Similarly, we find analytical expressions for the active flow velocity, pressure and vorticity related to the −1 defect, which are written in a compact form as The functions f a 3 , f a 5 , and f a ω giving the radial dependence of the velocity and vorticity are listed in the Appendix.These fields are also plotted in Fig. 3 showing that while the positive defect has closed streamlines also in an infinite system, the negative defect has an 8-fold symmetry of the vorticity, underlying a 4-fold symmetry of the polarity field around the defect.The streamlines do not close in an infinite system, but in confinement they might close due to boundary conditions or influence from other defects similarly to the vortices formed around ±1/2 defects in an active nematic 20,48 .It is important to note the distinctions of the analytical descriptions provided herein with a recent calculations of compressible flow fields inside circular confinements 37 .In the latter, the analyses is restricted to circular domains with radius smaller than the coherence length R 2 Kp/A, i.e the limit where relaxation/penetration length in the free energy is larger than the system size.Those results are therefore in the opposite limit of the calculation performed here.
In addition to the characteristic flow fields, our closed form analytic descriptions provide an insight into the isotropic stress patterns, i.e., half of the trace of the stress tensor, around full-integer topological defects.Such isotropic stresses have been shown to be a determining factor for the biological functionality of nematic defects, where concentration of compressive stress around +1/2 defects was shown to lead to cell death and extrusion 27 , while the tensile stresses around −1/2 defects has been shown to lead to spontaneous gap opening in epithelial cell layers 49 .Similarly, we find distinct isotropic stress patterns around positive and negative full-integer defects: while around a negative (−1) defect alternating regions of tension and compression appear in a 4-fold symmetric pattern (Fig. 3b), a strong augmentation of compressive stress is observed at the core of positive +1 defects (Fig. 2b).We conjecture that such a concentrated compressive stress at +1 defects could contribute to the activation of mechanosensitive signals in cell layer and potentially be linked to the recent observation of the cellular differentiation at +1 defect cores in cartilage cells 50 .

Vortical flow reversal around +1 defect
The results presented above provide closed-form analytical expressions for incompressible flows around full-integer topological defects in polar active matter.Interestingly, in the viscousdominated regime (Γ = 0), we find that the activity-induced flow around the +1 defect forms two counter-rotating vortices around the defect core, and thus the rotational flow reverses sign on a length scale comparable with the coherence length outside of the defect core.In this limit, the flow velocity becomes instead and its corresponding vorticity is thus where η = η/ξ 2 .We notice that the vorticity changes handedness for r = √ e ≈ 1.7 giving rise to two counter-rotating regions as shown in Fig. 4. Notice that the flow diverges at r → ∞ since there is no screening lengthscale.The mathematical reason for this flow reversal is that the logarithmic Greens function changes sign when r < 1, unlike the Bessel function K 0 kernel which appears when both friction and viscosity are present.Note that in the absence of friction, the flow equations have no natural lengthscale, thus the rescaling must be done with an appropriate lengthscale coming from the free energy, such as the coherence length which is consistent with the numerical observations of the flow reversal discussed next.Similar flow reversal is also expected in systems confined to circular domains that are small enough such that the energy scale can be treated as linearly dependent on the radial distance from the center of singularity χ(r) ∼ r in the entire domain 37 .For such systems the flow reversal happens inside the core region.

Comparison with numerical simulations
Here we discuss the comparison of the analytical predictions of active flow induced by ±1 defects with those obtained from direct numerical simulations of the polarization evolution Eq. (1-2) coupled with the Navier-Stokes equations with the additional passive stresses σp including the viscous stress as defined in Ref. 18 .In numerical simulations, we use a hybrid lattice-Boltzmann method, combining finite-difference for the evolution of polarity and the lattice-Boltzmann method for that of velocity.Using the same prescription as in Ref. 34  The colormap shows the magnitude of the velocity and pressure normalized by the maximum.Note that the pressure in the far-field diverges with the system size L. Since the pressure diverges at the centre the pressure inside the core r < 1 has been set to be constant.c) shows the normalised pressure, for an extensile system, at a fixed radius as a function of φ.Notice that for φ = ±π/4 + πn the pressure vanishes.
for density and viscosity, i.e. ρ = 40 and η = 3.6, we ensure that the Reynolds number in the simulations remains negligible (Re 1) 18,51 so that the dynamics of velocity virtually reduces to the incompressible Stokes equations Eq. (3-4) considered here.We fix the viscosity ratio to η/γ = 3.6, micro to macro length scale to ( Kp/A)/L = 2 × 10 −3 (assuring that the coherence length Kp/A is significantly smaller than the domain size L), and the flow alignment parameter to λ = 0.1.Dimensionless parameters in the simulations are defined similar to the theoretical parameters defined in Section 2.
Simulations were initialized with quiescent velocity field and noisy polar alignments close to the uniformly oriented state p = ex under periodic boundary conditions, on square domains of linear dimension L = 256.
We ran simulations with sufficiently large values of the activity parameter α0 for ±1 topological defects to form spontaneously 34 .We compile the flow profiles in the vicinity of each defect for a large sample of defects (ranging from 10 2 − 10 3 ) for each simulation.From these individual flow profiles, we compute the statistically-averaged flows generated by ±1 topological defects.However, because every defect had its own orientation, the flows had to be carefully reoriented upstream of the averaging procedure, as we detail below.
The orientation of −1 defects was defined from the surround-ing polarity field.Given the 2-fold symmetry of the −1 defect (see fig. 1d)), we search for the two principal axes along which the polarity vector points inwards and outwards, respectively, and align the local fields accordingly.Around +1 defects, both the polarity field and the velocity field are rotationally-symmetric (see Fig. 1a,b and Fig. 4).However, the vortical velocity around +1 defects may have different chirality: negative ("clockwise") when u θ < 0, or positive ("counterclockwise") when u θ > 0 as in Fig. 4. When simply summed up, these two chirality would cancel out, even though they have equivalent dynamic properties.To avoid this, we identify the flow chirality of each defect before averaging the fields.The flow chirality of each +1 defect was probed systematically by computing the sign of the circulation of u along small circular contours enclosing the defect center within the core of the defect, i.e. with the contour radius smaller than the length scale of the flow reversal √ e. Defects with clockwise chirality ( u θ < 0) are reversed (p θ ← −p θ , u θ ← −u θ ), while those with counterclockwise chirality are kept with the same orientation, so that all fields are consistently averaged.Note that no further distinction is made with regard to the polarity field, which therefore aggregates asters, spirals and vortices, since they lead to the same flow chirality.
The numerical results of the average polarity and flow profiles around ±1 defects are shown in Fig. 5 and Fig. 6.First, for both types of defects, the profile of the polarity field (right panels in Fig. 5 and Fig. 6) is consistent with the ideal vortex ansatz made by neglecting the backflow and the finite size of the core.The symmetries -of the average -are notably the same as in the passive case (see Fig. 1): positive defects are rotationally invariant -seemingly with a mixture of asters, spirals and vortices -, while negative defects feature a 2-fold symmetry.Nonetheless, since our model allows the polarity vector to shrink and vanish, the local magnitude of polarity constitutes an additional degree of freedom that yields qualitatively robust patterns around (−1)−defects, yet with a shift by a π/4−rotation between contractile ( ᾱ0 > 0) and extensile ( ᾱ0 < 0) systems, as can be seen from the colored fonts in fig.6( a, b and c, d ).These patterns become smaller as the net activity increases since the active lengthscale is inversely proportional to the activity l active ∼ K/α 0 (see fig. 6 a, c and b, d ).
We now discuss the flow fields generated by negative defects.As predicted analytically above, the negative full-integer defects feature 8-fold rotational symmetry, with regions where the velocity points alternatively inwards and outwards, interspersed with counter-rotating swirls.This structure is relatively robust, and strengthens as the level of activity increases.The flow structure is also very similar for contractile and extensile systems, up to a rotation by the same angle π/4− as above.In other words, equivalently, their profiles match, provided that one reverses the direction of the polarity vector according to the sign of the activity, which is consistent with the factor ᾱ0 in the expression of the field u we previously derived.Next, we discuss the velocity field generated by the positive defects.As expected from our theoretical considerations (fig.2, fig.4), the flow field exhibits a full-rotational symmetry.Since there was no friction, the dynamics was plainly viscousdominated, and we have indeed striking evidence of the flow reversal predicted by our calculations in this regime, for any level of activity: as in the previous fig.4, one can identify 2 inner and outer regions rotating contrariwise in fig. 5 a-d, that we will refer to as 'core' and 'shell' respectively.(In our captions, the core is rotating counterclockwise due to our averaging selection process).Taken together, the averaged velocity fields found from numerical simulation of polar active matter corroborate analytical predictions of the flow structures around full-integer defects and existence of counter-rotating vortices around positively charged defects.Furthermore, the simulation results were obtained from averaging the flow fields around full-integer defects within the active turbulence state and their close agreement with the analytical predictions around isolated defects give credit to the assumption of neglecting defect-defect interactions in the analytical calculations.

Defect pair interaction modulated by dipolar active forces
As seen in the above section the dipolar active force does not lead to any defect self-propulsion for isolated ±1 defects in an infinite domain.We now consider a pair of ±1 defects and show that the dipolar active force leads to an interaction between them through the active flow.For analytic tractability, we consider the frictiondominated regime Γ η and set the viscosity to zero, d → 0. The defect velocity is calculated using the methods from 38 (see 47 for details) and is given by for the +1 defect and −1 defect, respectively.Here Λ is a regularising lower cut-off corresponding to the finite defect core due to the divergent pressure at the +1 defect.R is the distance between the defects, ϕ is the angle of R = r − − r + .φ is the uniform background orientation field.This effective interaction induced by dipolar active forces decays inversely proportional with the distance between defects similar to the Coulomb-like force induced by the phase gradients alone.However, these interactions are anisotropic since they depend on the orientation of the background polarization as well as the orientation of the defect pair.

The effect of polar active forces
Next, we also include the polar active force in the main flow equations Eq. ( 5) and study its contribution to the flow velocity up and pressure Pp.In dimensionless units, the polar active force is where αp = αpτ /(ξΓ) is a rescaled parameter that measures the strength of polar active forces relative the frictional drag.We consider αp > 0 corresponding to polar particles moving in the direction of their head.The corresponding flow velocity and pressure fields can now be written compactly as (see details in SM 47 ) for the positive and the negative defect respectively.The velocity profiles are plotted in the Fig. 7 while the functions f 1 and f 3 are plotted in Fig. 8 together with their asymptotic values.It is straightforward to show that both f 1 and f 3 tend to zero at the defect origin, thus the polar active forces do not contribute to defect self-induced motility.It is important to note that the flow induced by polar active forces around the +1 defect is proportional to sin(φ), instead of the sin(2φ) dependence that was found for   Fig. 8 Plot of the functions (blue, whole lines) f p 1 (a) and f p 3 (b) against their asymptotic limits (black, doted lines) as given in eq. ( 32) and (33).Note that the functions asymptotically approach to non-zero values in the far-field.
the flow induced by dipolar active forces (Eq.16).This means that while the polar active forces do not produce any flow around an ideal aster (φ = 0), there will be a polar activity-induced flow around an ideal vortex (φ = π/2).Furthermore, Eq. 29 shows that the pressure is proportional to the radial distance ∼ r, indicating that, in the absence of any far-field screenings from other defects, the pressure can become quite large.This is related to the flow induced by polar active forces being constant in magnitude in the fare-field, as illustrated in Fig. 8, or by explicitly writing the asymptotic limits for the negative and positive defect respectively.As we will show explicitly for a pair of defects in the next section, Eq. 32-33, together with Eq. ( 13), indicate that there are strong interactions between defects regardless of the distance between them.Having found closed form formulas for the velocity field around the defects, the far-field vorticity induced by the positive defect can be calculated as and, for the negative defect, this is Interestingly the vorticity of the positive defect is independent on the polar angle θ, which is also true in the near-field, because we can write the velocity field in the form u = r ⊥ f (r).As such, close to the defect center the vorticity of the negative defect vanishes at the center, while the positive defects has a finite vorticity indicating that the positive full-integer defect will be endowed with a spin.

Defect-pair interaction in friction-dominated system
Taking again the analytically tractable limit of zero viscosity, we can find the velocity induced by polar active forces for a pair of oppositely-charged defects following the approach from Sec. 3.3.The flow velocities in the centre of the ±1 defects reduce to where R = r − − r + is the separation vector between the negative defect and the positive one, and p 0 = cos(φ) ex + sin(φ) e − is the uniform background polarity.Interestingly, the polar active forces induce non-reciprocal and non-local mutual interactions that depend only on the orientation of R relative to φ, independent of the separation distance R between the defects.This suggests a truly long-ranged interaction between oppositely-charged defects in the presence of polar active forces.Furthermore, the −1 defect tends to move towards/away from +1 depending on the orientation angle of R with respect to p.Both defects move perpendicular to R but at different rates inducing the pair to rotate.We can see this behavior more clearly by looking at the defect pair velocity under polar active forces alone, where the first term ( R) is an attraction/repulsion force between the defects, while the second term (⊥ R) rotates the defect pair.The rotation is zero, when the defect pair aligns with the background polarization in either directions ( R p 0 ).We show in SM. 47 that the defect pair rotates until annihilation, unless R is initially parallel to p 0 .This shows that the polar active force leads to pair annihilation of oppositely-charged defects to promote large-scale polar order.
To check this interesting dynamics numerically, we compare the pair trajectory determined by Eq. ( 13) with that predicted by the full hydrodynamic model from integrating Eqs.(1) (3) and ( 4) with the same initial configuration.The initial uniform polarisation p = ey is seeded with a defect pair in the x-direction, i.e. r − = −20 ey and r + = 20 ex, such that R ⊥ p 0 at t = 0.The model parameters are set to αp = 0.5, α0 = 0, λ = 0 and ζ 2 = 0.01.The hydrodynamic equations are solved with periodic boundary conditions in a domain 256×256 with a spatial discretisation ∆ = 0.5, using spectral methods and an exponential time differentiation scheme 52 .Under the polar active force, the separation vector R rotates relative to p 0 , and this changes the shape of the +1 defect from an initial vortex to an aster at the annihilation time.This effect is lost in the absence of polar active forces, i.e.at αp = 0 (see animations of the defect pair annihilation in 47 ).
The two trajectories are shown in Fig. 9a) and agree very well for large R = | R| in consistency with the pointwise approximation.It also shows that polar active forces increase the annihilation rate.This is further evidenced in Fig. 9 b) where we compare the evolution of the separation distance R(t) when the two defects interact through the Coulomb-like forces with or without the presence of polar active forces.At αp = 0, R(t) decreases in the far-field as √ t 0 − t (t 0 being the annihilation time) due to the 1/R interaction forces.However, for αp = 0, the annihilation timescale is greatly reduced and the deviations from the √ t 0 − t behavior arise due the non-local and non-reciprocal attraction force.The rotation rate is also contributing to aligning the defects to increase the attraction between the defects.

Discussion/conclusion
In summary, we present theoretical derivations of incompressible flow fields and dynamics of ±1 defects under dipolar and polar active forces.These active forces do not endow the full-integer defects with any self-propulsion, as expected from the symmetry of the defects.However, the +1 defect acquires a non-zero active torque due to both polar and dipolar active forces.We show that the strength of this spin is dependent on whether the defect is a vortex, aster or spiral.In the absence of hydrodynamic screening due to friction, the vortical active flow around a +1 vortex changes sign on a length scale set by the coherence length as observed in numerical simulations and predicted analytically.For both defects, the dipolar active force contribution to the flow field vanishes with 1/r in the defect far-field, whereas the polar active force contribution approach a constant value.Remarkably, polar active forces mediate mutual interaction between oppositely-charged defect pairs in a manner that renders defectdefect interactions independent of the distance between the defect pair.We have shown that, under these polar active forces, a pair of oppositely-charged defects rotates to align with the background polarization field, while the negative defect chases the positive one until annihilation.The rate of annihilation is greatly enhanced by polar active forces, and this is the main underlying mechanism for the suppression of defect-laden active turbulence in polar active matter as reported numerically in Ref. 34 .

Appendix:
Here we write down the radially-dependent functions associated with the flow fields induced by polar and dipolar active forces.Derivation details are found in the SM. 47.The radial functions for the dipolar active flow velocity are: The radial function for the corresponding vorticity field f a ω (r) = 1 r4 −48 + 4r 2 + r2 (24 + r2 )K 0 (r) + 8r(6 + r2 )K 1 (r) . ( The radial functions for the polar active flow velocity are: and Where ψ (0) is the digamma function.The coefficients of the infinite power series .(48)   Note that when summing over n the coefficients κ 4 (n, k) tends to zero for any k.The infinite series are convergent, but slowly.

Fig. 1
Fig. 1 Polarity field around the ±1 topological defects.a) and b) show a +1 defect, while c) and d) illustrate a −1 defect.The drawings a) and c) show the p plotted on a circle with the angle φ marked.This illustrates why φ can be set to zero by a change of basis for the −1 defect and not for the +1 defect.b) and d) show the polarity fields around a +1 defect with φ = π/3 (spiral) and a −1 defect with φ = 0 respectively (Changing φ rotates the negative defect).

Fig. 2
Fig.2Incompressible active flow streamlines (a) and corresponding pressure field (b) around an isolated +1 vortex with sin(2φ) > 0 in an extensile ( α0 < 0) system.The colormap shows the magnitude of the velocity and pressure normalized by the maximum.Note that the pressure in the far-field diverges with the system size L. Since the pressure diverges at the centre the pressure inside the core r < 1 has been set to be constant.c) shows the normalised pressure, for an extensile system, at a fixed radius as a function of φ.Notice that for φ = ±π/4 + πn the pressure vanishes.

Fig. 3
Fig. 3 Incompressible active flow streamlines (a) and corresponding pressure field (b) around an isolated −1 vortex with φ = 0 in an extensile ( α0 < 0) system.The colormaps show the normalised magnitude of velocity (a) and pressure (b).

Fig. 4
Fig. 4 (a) Streamlines of the theoretical active flow around an isolated +1 vortex for sin(2φ) > 0 for an extensive system and Γ = 0.The colormap shows the flow magnitude.The change in color of the streamlines from white to black highlight the flow reversal r ≈ 1.7.This is more apparent in (b) which is a cross section of the azimuthal velocity component along the positive x axis.The velocity changes sign when it crosses the zero line (doted black).

1 Fig. 5
Fig. 5 Numerical average polarity and velocity fields around (+1)−defects, for various activities.p (resp.u) field is plotted on the right (resp.left) of each panel.Averages were computed on 100 defects for α0 = ±0.4,and on 1000 defects for α0 = ±1.The polarity field consists each time of a statistical combination of asters, spirals and vortices.

− 1 Fig. 6
Fig. 6 Numerical average polarity and velocity fields around (−1)−defects, for various activities.p (resp.u) field is plotted on the right (resp.left) of each panel.Averages were computed on 100 defects for α0 = ±0.4,and on 1000 defects for α0 = ±1.Abnormal points (max.8, on the edge of the averaging window) were set down to 0 to retain only significant flow patterns.

Fig. 7
Fig. 7 Velocity streamlines due to polar active forces for a) the positive defect (vortex flow) and b) the negative defect (saddle-point flow).Colrepresents the magnitude of the velocity field normalized by its maximum value.
) J o u r n a l N a me , [ y e a r ] , [ v o l .] , 1-14 | 9

Fig. 9 a
Fig. 9 a) Phase portraits of the defect pair trajectory obtained from (blue) the kinematic law for point defects (PD) from Eq. (13), and (red) by direct simulation of the full hydrodynamic model (HM) starting from the same initial condition.The positive defect is marked by +, the negative by −. αp = 0.5, α0 = 0, λ = 0 and ζ 2 = 0.01.b) Distance R between a defect pair as a function of time for Coulomb-like interactions only (blue curve) and in the presence of polar active forces (orange curve).The black doted lines are fitting curves with the R(t) = √ At + B found by linear regression on R 2 .