††thanks: These authors contributed equally to this work.††thanks: These
authors contributed equally to this work.††thanks: corresponding author:
<EMAIL_ADDRESS>
# Quadrupolar active stress induces exotic phases of defect motion in active
nematics

Salik A. Sultan  The Niels Bohr Institute, University of Copenhagen,
Copenhagen, Denmark  Mehrana R. Nejad  Rudolf Peierls Centre for Theoretical
Physics, University of Oxford, UK  Amin Doostmohammadi  The Niels Bohr
Institute, University of Copenhagen, Copenhagen, Denmark

###### Abstract

A wide range of living and artificial active matter exists in close contact
with substrates and under strong confinement, where in addition to dipolar
active stresses, quadrupolar active stresses can become important. Here, we
numerically investigate the impact of quadrupolar non-equilibrium forces on
the emergent patterns of self-organisation in non-momentum conserving active
nematics. Our results reveal that beyond having stabilising effects, the
quadrupolar active forces can induce various modes of topological defect
motion in active nematics. In particular, we find the emergence of both polar
and nematic ordering of the defects, as well as new phases of self-
organisation that comprise topological defect chains and topological defect
asters. The results contribute to further understanding of emergent patterns
of collective motion and non-equilibrium self-organisation in active matter.

##  I Introduction

Active materials are out of equilibrium systems in which energy is
continuously consumed and transferred to mechanical work by constituent
elements of the matter Sanchez  _et al._ (2012); Weirich  _et al._ (2019);
Marchetti  _et al._ (2013); Needleman and Dogic (2017); Jülicher  _et al._
(2018). Examples, at various length scales, include flocks of birds Toner and
Tu (1998), schools of fish, networks of actin and microtubules Bate  _et al._
(2019) and cellular tissues Duclos  _et al._ (2018); Kawaguchi  _et al._
(2017). In contrast to most non-equilibrium systems, where an external driving
force is applied to the system, the defining attribute of active materials is
that the energy is input locally at the level of individual components of a
structure, leading the system into a dynamic state. The local energy injection
results in the generation of active stresses by the particles that in turn
drive distinct patterns of collective motion and self-organisation of active
materials. Examples range from large-scale coordinate motion of cells Poujade
_et al._ (2007) and bacteria Meacock  _et al._ (2020), to coherent flows of
subcellular filaments under geometric confinement Wu  _et al._ (2017).

The active stresses, therefore, are often considered to induce destabilising
effects, putting an otherwise quiescent matter into motion Aditi Simha and
Ramaswamy (2002); Voituriez  _et al._ (2005); Thampi  _et al._ (2015);
Giordano  _et al._ (2021). The destabilising effects of active stresses can
even lead to chaotic flows of active particles in a state termed active
turbulence Giomi  _et al._ (2011); Thampi and Yeomans (2016). However, because
of their widespread implications in biological processes and in the design of
non-equilibrium materials, there is a growing interest in control and
stabilisation of active materials Needleman and Dogic (2017); Doostmohammadi
_et al._ (2018); Brusatin  _et al._ (2018). In this venue, several theoretical
and experimental mechanisms have been proposed to tame, the otherwise chaotic
motion of active particles Doostmohammadi and Yeomans (2019); Shankar  _et
al._ (2021). Particularly successful approaches use topological Keber  _et
al._ (2014); Metselaar  _et al._ (2019); Carenza  _et al._ (2020) and
geometric constraints Hardoüin  _et al._ (2019) to streamline flows of active
particles. Mechanistically, such constraints induce hydrodynamic screening
effects that allows for stabilisation of active flows. For example, placing
bacterial suspensions or microtubule-motor protein mixtures under confinement
is shown to result in a crossover from chaotic flows to vortex-lattices and
coherent streams as the confinement size decreases Wioland  _et al._ (2016);
Hardoüin  _et al._ (2019). Similarly, placing active particles in contact with
substrates results in friction-induced screening length that could stabilise
chaotic flows into vortex-lattices Doostmohammadi  _et al._ (2016); Thijssen
_et al._ (2020).

Not only do these screening effects extract momentum from the system, they
also affect the dominant forms of active stresses that are exerted by
particles on their surroundings. In the absence of screening effects, it is
well-established that the dominant stress contribution from a mirco-scale
active particle is a stresslet Mak  _et al._ (2016); Ramaswamy (2010); Prost
_et al._ (2015). However, when an active particle is in contact with a
substrate or under extreme confinement, the dipolar stress is screened and a
higher order quadrupolar contribution also becomes important Mathijssen  _et
al._ (2016). In this vein, a recent study has introduced the coarse-grained
continuum representation of quadrupolar forces as an additional non-
equilibrium active force that is generated by dense suspension of active
particles Maitra  _et al._ (2018). Interestingly, using linear stability
analyses, it is shown that such quadrupolar force can induce a stabilising
effect on the active system, counteracting the hydrodynamic instability that
is induced by the dipolar forces. Notwithstanding this important contribution,
the non-linear impacts of such additional non-equilibrium force on the pattern
formation and dynamic organisation of active systems remain unexplored. Going
beyond the linear instability, here we numerically examine the full non-linear
effects of the additional quadrupolar force on the dynamics of collective
pattern formation in a model of active nematics. We focus particularly on the
dynamics of the topological defects, singularities in the orientation field,
that are increasingly identified in various realisations of active matter and
are shown to have biological functions Saw  _et al._ (2017); Kawaguchi  _et
al._ (2017); Copenhagen  _et al._ (2020); Meacock  _et al._ (2020); Maroudas-
Sacks  _et al._ (2020); Zhang  _et al._ (2021); Blanch-Mercader  _et al._
(2021); Doostmohammadi and Ladoux (2021). The results uncover exotic phases of
defect organisation that comprise topological defect chains, asters of
topological defects, as well as polar and nematic ordering of the motile
topological defects.

Figure 1: Phase space of the system for different values of the tumbling
parameter $\chi$ and the activity ratio $\mathcal{Z}=\zeta_{2}/\zeta_{1}$. The
nematic order can become stable for positive values of the force quadrupole
term. The yellow region in the phase diagram shows the stable region found
analytically in Ref. Maitra  _et al._ (2018). Snapshots showing the
distinguishable characterization of defect chains, polar defect ordering,
defect asters, polar defect flocks. The underlying color map shows the
magnitude of the nematic order, and $+1/2$ ($-1/2$) defects are represented in
green (blue).

The paper is organised as follows: in section II we introduce our model and
discuss the details of the governing equations, active forces, and provide the
list of simulation parameters. In section III, we introduce different phases
that are revealed through numerical simulations of the system and provide a
phase diagram that represents the emergent phases as a function of dipolar and
quadrupolar force strengths and the flow-aligning parameter. In section IV we
represent a quantitative characterisation of observed phases and analyse
emergent topological defect organisations. Finally, concluding remarks and
future directions are discussed.

##  II Methods

In order to understand mechanisms behind diverse behaviours of active
materials, it is often useful to divide active systems into two general
categories: wet and dry active matter. In wet active matter, flows created by
the particles mediate long-range hydrodynamic interactions in the system Koch
and Subramanian (2011) such as in unconfined suspensions of filaments and
motor proteins Needleman and Dogic (2017); Kumar  _et al._ (2018), bacterial
suspensions in low concentrations Turiv  _et al._ (2020), and collection of
artificial active colloids suspended within a fluid medium Bäuerle  _et al._
(2020). On the other hand, the collective behaviour of dry active materials is
governed by direct interactions, such as collision between active particles
Scholz  _et al._ (2018); Chaté (2020). Shaken granular matter Junot  _et al._
(2017); Kumar  _et al._ (2014) and bacteria moving on a substrate Meacock  _et
al._ (2020) are examples of dry active matter. Descriptions based on dry
active matter theories have been helpful in investigating the different
patterns in clusters of bacteria, which live on dry surfaces in tight spaces
Peruani  _et al._ (2012); Großmann  _et al._ (2016). Since in this work our
focus is on the behavior of active model system in the presence of strong
hydrodynamic screening, we employ equations of dry active nematics
Doostmohammadi  _et al._ (2016) with the additional quadrupolar active forces.
The active nematic model is chosen since it has proven successful in
describing a wide variety of active materials that constitute elongated
building blocks such as rod-shaped bacteria Nagai (2018), subcellular
filaments Alvarado  _et al._ (2014), and spindle-shaped cells Duclos  _et al._
(2017), as well as deformable cells in which orientational ordering is an
emergent feature Mueller  _et al._ (2019); Saw  _et al._ (2017); Blanch-
Mercader  _et al._ (2018).

To describe the orientational order, we study the evolution of the nematic
order parameter $\textbf{Q}=2q(\hat{n}\hat{n}-\frac{\textbf{I}}{2})$, where
$q$ and $\hat{n}$ show the magnitude and orientation of the order,
respectively, and I represents the identity tensor.
The dynamics of the Q-tensor is coupled to the velocity of the surrounding
fluid $\bm{u}$, and reads:

$\displaystyle\frac{DQ_{ij}}{Dt}$ $\displaystyle=\chi E_{ij}+\Gamma H_{ij},$
(1) $\displaystyle fu_{i}$ $\displaystyle=-\zeta_{1}\
\partial_{j}Q_{ij}-\zeta_{2}Q_{ij}\partial_{k}Q_{kj},$ (2)

where we have introduced the convected co-rotational time derivative
$DQ_{ij}/Dt=(\partial_{t}+u_{k}\partial_{k})\
Q_{ij}-\omega_{ik}Q_{kj}+Q_{ik}\omega_{kj}$, in which
$\omega_{ij}=(\partial_{j}u_{i}-\partial_{i}u_{j})/2$ is the vorticity tensor.
In Eq. 1, $\chi$ is the tumbling parameter and
$E_{ij}=(\partial_{i}u_{j}+\partial_{j}u_{i})/2$ is the strain rate tensor.
$\Gamma$ shows the rotational diffusivity which along with $H_{ij}$, the
molecular field, describes the relaxation of the Q-tensor towards the minimum
of a free energy described by:

$\mathcal{F}=\frac{A}{2}(q^{2}-\frac{1}{2}\textbf{Q}:\textbf{Q})^{2}+\frac{K}{2}|\nabla\textbf{Q}|^{2}+g(\nabla\cdot[\nabla\cdot\textbf{Q}])^{2}.$
(3)

The free energy includes a nematic alignment term (with coefficient $A$) and
an elastic term which penalises gradients in the Q tensor. Here, we use a
single elastic constant approximation and show the elastic constant by $K$. In
line with previous studies of dry active nematic Oza and Dunkel (2016); Putzig
_et al._ (2016), we also include a regularization term with coefficient $g$
which provides stability to our numerical solutions at small length scales.

In Eq. 2, $f$ is friction coefficient, and $\zeta_{1}$ and $\zeta_{2}$ show
activity coefficient related to a force dipole and a force quadrupole,
respectively. Two well-known classes of active suspensions, namely extensile
and contractile systems, are described by $\zeta_{1}>0$ and $\zeta_{1}<0$,
respectively. Extensile systems comprise elements that pump the surrounding
fluid in an outwards direction along their symmetry axis $\hat{n}$ and inwards
in the direction perpendicular to $\hat{n}$. Contrarily, in a contractile
system, elements pump fluid inwards along $\hat{n}$ and outwards in a
direction perpendicular to that. As a result of the flows created in active
suspension in the bulk, an extensile system is stable under splay deformation
and unstable under bend deformation Aditi Simha and Ramaswamy (2002).
Contrarily, a contractile system is stable under bend deformations, but
unstable as a result of splay deformations Thampi and Yeomans (2016). In
addition to this well-established dipolar activity, the force distribution of
the higher angular symmetry, the force quadrupole term with coefficient
$\zeta_{2}$ in Eq. 2, becomes important in systems in contact with a
substrate, where the momentum is not conserved Maitra  _et al._ (2018).

We study the role of the quadrupole active force in the dynamics of the
system, by numerically solving the coupled equations of the nematic tensor Q
and the velocity field $\bm{u}$. In this regard, we fix the value of the
dipolar active force to $\zeta_{1}=0.2$, and study the role of the tumbling
parameter $\chi$ as well as the strength of the active force quadrupole
$\zeta_{2}$ in the dynamics. As such, we vary the strength of the active
quadrupole term in the range $-0.3\leq\zeta_{2}\leq 0.25$, effectively
exploring the interplay between the dimensionless ratio
$\mathcal{Z}=\zeta_{2}/\zeta_{1}$ and the tumbling parameter $\chi$, which
determines whether the system is in a flow tumbling or flow-aligned regime
Ternet  _et al._ (1999). In a flow tumbling regime, the director tumbles under
shear flow while in a flow aligning regime it aligns with an angle (Leslie
angle) with the direction of the shear Beris and Edwards (1994). The system
size of $512\times 512$ grid points is simulated with periodic boundary
conditions on all sides. Parameters for the rotational diffusivity,
regularisation, and friction are fixed at $\Gamma=0.05$, $g=0.1$, and
$\textit{f}=10$. In addition, the elasticity, and the bulk nematic alignment
coefficient are set to $\text{K}=0.15$ and $\text{A}=0.5$, respectively.

##  III Results

Figure 2: a) Root mean square velocity $v_{\text{rms}}$ averaged over space
and time in steady state. Regardless of flow-aligning or tumbling state
$v_{\text{rms}}$ decreases as $\mathcal{Z}$ increases. By proxy, it is also
apparent that increasing the quadrupolar active term, lowers the velocity of
the system. b) Magnitude of the order defined as $\langle
Q_{xx}^{2}+Q_{yx}^{2}\rangle$ averaged over space and time in steady state.
The average order follows the inverse trend compared to the rms-velocity, and
increases by increasing $\mathcal{Z}$. Figure 3: Average deformation defined
as $\langle\Delta
n\rangle=\langle(\partial_{i}\partial_{j}Q_{ij})^{2}\rangle$. The average is
taken over time and space in steady state The fall in deformation correlates
with observed stabilisation as a consequence of increasing $\mathcal{Z}$.

###  III.1 Stability diagram

We begin by investigating the emergent patterns for varying values of the
quadrupolar to dipolar activity ratio $\mathcal{Z}=\zeta_{2}/\zeta_{1}$ and
the flow-aligning parameter $\chi$. The results are summarised in a stability
diagram in Fig. 1, where distinct patterns of motion are represented.

Fig. 1 shows that a positive value of the quadrupolar activity can stabilise
the nematic phase and the required quadrupolar force for stability increases
with increasing the value of the flow-aligning parameter. This is in agreement
with the analytical result of Ref. Maitra  _et al._ (2018): in Fig. 1 we show
the stable nematic phase predicted analytically in yellow background color.

At the border between stable nematic phase and unstable phases, flows induced
by dipolar activity are just strong enough to trigger nematic deformations but
not strong enough to nucleate topological defects and instead we observe the
formation of _extensile walls_. These extensile walls have also been observed
previously Nejad  _et al._ (2021) in wet active nematics in the presence of
friction, and are a consequence of the competition between active forces to
create defects and elastic forces to maintain the nematic phase. Going beyond
the stability boundary, topological defects are nucleated in the system and we
observe the emergence of several dynamic phases of defect organisation as a
result of the competition between dipolar and quadrupolar active forces.

To quantitatively distinguish between the stable phase (Fig 1, yellow region),
where we either observe a nematic phase or extensile walls, and an unstable
phase, where topological defects form and patterns emerge, we calculate the
root-mean-square velocity averaged over both space and time after the
statistical steady-state is reached. As evident from Fig. 2a active flows
decrease with increasing the quadrupolar force coefficient $\zeta_{2}$, which
is consistent with the stabilising effect of the quadrupolar term.

The development of active turbulence stems from the onset of instabilities in
the director field. Parallel walls consisting of lines of bend deformation
separate the nematic regions. The decay of these walls leads to the formation
of topological defects Thampi  _et al._ (2014). In order to quantify the
nematic order and distinguish the nematic phase, we measure the average
magnitude of the nematic order (defined as $S_{0}=\langle
Q_{xx}^{2}+Q_{yx}^{2}\rangle$). This quantity measures the global nematic
order and thus indicates the overall stability of the nematic phase. A
decrease in the value of the nematic order indicates the formation and
presence of more defects. The average order is shown in Fig. 2b for different
values of the tumbling parameter. These graphs show the same trend as the
average mean root square velocity graphs, indicating that flows are mainly
created in places with large drop in the magnitude of the nematic order. Both
graphs show a transition from the disordered phase to ordered phase by
increasing the quadrupole force.

Figure 4: Time evolution of dynamics for (a) polar defect flock and (b) defect
chains. In the polar defect flocking phase, the system forms motile defects
that flock together in anti-parallel directions as indicated by the arrows.
Defects in the chain phase self-organise to form (c) long chains of robust
$-1/2$ defects as well as (d) chains of dynamic oscillating $+1/2$ defects as
the system relaxes towards a steady state.

To further distinguish between the ordered and disordered phases of defect
configurations we measure the average deformation, defined as $\langle\Delta
n\rangle=\langle(\partial_{i}\partial_{j}Q_{ij})^{2}\rangle$, for different
values of the tumbling parameter and as a function of activity ratio
$\mathcal{Z}$. A large value of deformation follows a distortion in Q-tensor.
Note that this quantity is different from the average order $S_{0}$ as it
accounts for the distortions in the director and the deviation from the global
nematic ordered phase. $\langle\Delta n\rangle$ is maximum in the presence of
walls or in the active turbulence phase, as in these phases the global
distortions in the director is maximum. The average distortion is presented in
Fig. 3 for varying values of the tumbling parameter. Note that ordered phases
of defect organisation do not correspond to a peak in the distortion graph as
in these phases the distortions in the director are local (in places of
defects) and not global.

Next, we describe the characteristics of defect self-organisation as this new
non-linear regime is the focus of the current study.

###  III.2 Dynamic phases

We now turn to the unstable part of the phase diagram, where for a given value
of the flow-aligning parameter, the instability driven by dipolar activity is
strong enough to lead to deformation of the director field and nucleation of
topological defects, while the presence of the quadrupolar activity acts to
restructure the spatial defect organisation. In this regime, we found a
variety of phases depending on the value of the tumbling parameter:

Defects flocking: In the unstable phase, but close to the boundary of the
stable phase, we observe a phase in which different groups of self-propelled
$+1/2$ defects move together forming flocks of motile topological defects
(Fig. 1 and Supplementary Movie $1$ in sup ). Fig. 5a shows the angular
distribution of $+1/2$ defects within this phase. Defects are identified using
the diffusive charge density Blow  _et al._ (2014)

$s=\left(\frac{1}{2\pi}\right)\times\left(\left[\frac{\partial
Q_{xx}}{\partial x}\frac{\partial Q_{xy}}{\partial
y}\right]-\left[\frac{\partial Q_{xx}}{\partial y}\frac{\partial
Q_{xy}}{\partial x}\right]\right),$

which takes values of $\pm 1/2$ at the defect cores and is zero otherwise. The
angle $\theta$ is calculated using
$\nabla\cdot\mathbf{Q}=(\cos\theta,\sin\theta)$, where $\nabla\cdot\mathbf{Q}$
gives the direction of the defects self-propulsion Pearce  _et al._ (2021)
that for extensile systems ($\zeta_{1}>0$) is from tail to head of the $+1/2$
defect. Interestingly, in this phase, the defects in each flock move in the
same direction and show polar order, but different defect flocks can migrate
in anti-parallel directions within the system as seen from the snapshots in
Fig. 4a. To show this behavior more clearly, we calculate both polar and
nematic order parameter for the alignment of the $+1/2$ topological defects.
The polar ordering is defined by

$P_{\text{+1/2}}=\sqrt{P_{x}^{2}+P_{y}^{2}},$ (4)

where

$[P_{x},P_{y}]=\left[\sum_{i}^{N}\frac{m_{xi}}{N},\sum_{i}^{N}\frac{m_{yi}}{N}\right],$
(5)

and $m$ is defined by the orientational angle $\theta_{i}$ of each $+1/2$
defect as $m_{xi}=\cos\theta_{i}$ and $m_{yi}=\sin\theta_{i}$. Similarly, the
nematic order of the defect alignment is characterised by the higher order
multipole of the alignment angle De Gennes and Prost (1993)

$S_{\text{+1/2}}=\frac{1}{2}\langle 3\cos^{2}\theta-1\rangle.$ (6)

For a perfect polar alignment of the motile defects $P_{\text{+1/2}}=1.0$,
while for the perfect nematic alignment $S_{\text{+1/2}}=1.0$. For the angular
distribution represented in Fig. 5a that corresponds to the defects flocking
state, we find $P_{\text{+1/2}}=0.24$ and $S_{\text{+1/2}}=0.8$, indicating
that the orientational organisation of the $+1/2$ defects in the defect flocks
state shows a dominant nematic alignment.

Polar defect ordering: Close to the defect flocking state, reducing the flow-
aligning parameter together with lower quadrupolar activity leads to the
global polar ordering of the motile $+1/2$ defects (Fig. 1 and Supplementary
Movie $2$ in sup ). This is reminiscent of the polar defect ordering predicted
by the hydrodynamic theory of active defects in the over-damped limit and is
due to the active aligning torque on defects in a charge neutral system
Shankar and Marchetti (2019). The polar ordering state is best represented in
the angular distribution plot showing the sharp peak in the orientational
alignment of the defects (Fig. 5b) and corresponds to the quantitative value
of the global polar order $P_{\text{+1/2}}=0.5$.

Defect asters and active turbulence: At the zero value of the flow-aligning
parameter, where the director only responds to the rotational part of the flow
gradient, and for negligible values of the quadrupolar activity coefficient,
we find a different state characterised by aster-like organisation of groups
of $+1/2$ defects (Fig. 1 and Supplementary Movie $3$ in sup ). Increasing the
strength of the quadrupolar activity towards negative values results in the
breakup of the aster-like structures and the establishment of the active
turbulence characterised by the chaotic motion of the topological defects.
Within both phases the orientational organisation of the defects is isotropic
(Fig. 5c,d and Supplementary Movie $4$ in sup ).

Figure 5: Angular distribution of $+1/2$ defects for different phases: a)
Defects flocking, b) Polar defect order, c) Active turbulence d) Defect
asters, e) Defect chain and f) Polar defect flock coexisting with the defect
chain. The different phases are distinguishable both by the direction of the
$+1/2$ defects as well as the number of defects that form in a given phase.
The colourbar corresponds to the number of $+1/2$ defects orientated at angle
$\theta$

Defect chains: Interestingly, for negative moderate quadrupolar activities, we
found a phase in which $\pm 1/2$ defects self-organise into large chain-like
structures. The chains of $+1/2$ and $-1/2$ defects show different
characteristics: while $-1/2$ defect chains are stable and do not restructure
after they are formed, the $+1/2$ defect chains are dynamic and within the
chains pairs of $\pm 1/2$ defects nucleate, annihilate, move, and oscillate
around each other (see Fig. 4b and Supplementary Movie $5$ in sup ) in a
fashion resembling local defect dancing that has been reported for confined
active nematics Doostmohammadi and Yeomans (2019). The large-scale defect
chains could be reminiscent of the propagating soliton-like clusters that have
been recently discovered for confined active nematics with chiral anchoring Li
_et al._ (2021), though here they are emergent properties of the system, where
neither confinement nor chiral anchoring are present. In this vein, the impact
of the friction on hydrodynamic screening and quadrupolar activity on breaking
the angular momentum conservation could resemble the confinements and chiral
anchoring effects, respectively, that have been shown to induce soliton-like
defect clusters.

In order to better understand the mechanism of chain formation we analytically
calculate the force induced by dipolar and quadrupolar activities around
isolated $\pm 1/2$ topological defects. The Q tensor around a topological
defect with charge $k$ reads:

$\textbf{Q}=\begin{pmatrix}\cos(2k\phi)&\sin(2k\phi)\\\
-\sin(2k\phi)&\cos(2k\phi)\end{pmatrix},$ (7)

where the defect symmetry axis is along the $x$-axis (as represented in Fig.
6), and $\phi$ is the polar angle (measured counter-clockwise with respect to
the $x$-axis). Using this definition for the Q tensor the active force around
a defect reads:

$\textbf{f}^{a}=\frac{qk}{r}\begin{pmatrix}\zeta_{2}\cos{\phi}-\zeta_{1}\cos{\phi(1-2k)}\\\
\zeta_{2}\sin{\phi+\zeta_{1}\sin{\phi(1-2k)}}\end{pmatrix}.$ (8)

From the symmetry of the quadrupolar active force it is expected that such a
force does not have any impact on the self-propulsion speed of the defects.
Interestingly, however, the quadrupolar activity induces distinct diverging
and converging forces around $+1/2$ and $-1/2$ topological defects,
respectively (Fig. 6). The converging force around an isolated $-1/2$ defect
can explain the attraction of $-1/2$ defect pairs and the formation of chains
of $-1/2$ defects (Fig. 4b,c). On the other hand, the $+1/2$ defects that are
in small distance from the $-1/2$ chain annihilate, but the $+1/2$ defects far
from the $-1/2$ chain can form a dynamic chain. The orientation of the defects
in a $+1/2$ chain, agrees with the active torque between $+1/2$ defects
introduced in Ref. Shankar  _et al._ (2018) and suggests that accumulation of
$-1/2$ defects in chains allows $+1/2$ defects to experience an aligning
active torque purely caused by the other $+1/2$ defects.

Finally, deep in the unstable phase we found a phase in which the polar defect
flock phase coexists with the defect chain phase. In both pure defect chain
phase and the coexistence phase, characterisation of the angular defect
orientation reveals a dominant nematic alignment of the $+1/2$ defects (Fig.
5e, f).

Figure 6: Force distribution around $\pm 1/2$ topological defects due to the
dipolar and quadrupolar active forces, showing distinct diverging and
converging force patterns due to the quadrupolar force around $+1/2$ and
$-1/2$ defects, respectively. Here, blue solid lines indicate the director
field of the defects, while yellow and red arrows illustrate force induced by
dipolar and quadrupolar activities, respectively. Here, we have used
$\mathcal{Z}=-0.75$.

##  IV Conclusion

We have numerically studied the dynamics of active nematic systems in contact
with a substrate. In this setup, a new active force contribution, that is
absent in a bulk momentum conserving system, plays a role in the dynamics of
the system Maitra  _et al._ (2018). We confirmed that for one sign of this
active force, the nematic state becomes stable in agreement with the
analytical result of Ref. Maitra  _et al._ (2018). Remarkably, going beyond
linear stability, we also showed that in the unstable regime where the new
active force cannot recover the nematic phase, it introduces new patterns of
topological defect organisation in the system. In particular, we found the
formation of stable chains of $-1/2$ topological defects, while the
corresponding $+1/2$ defects form elongated clusters and oscillate around a
center in pairs. By varying the tumbling parameter, we additionally found the
emergence of polar defect flocks, in which system forms different flocks of
$+1/2$ defects. Defects in each flock represent polar order but different
flocks move in opposite direction and show nematic order. Furthermore, states
showing pure polar flocking of the $+1/2$ defects and aster-like dynamic
defect configurations, in which $+1/2$ defects point radially towards a
center, were uncovered for varying quadrupolar-to-dipolar active force ratios
and flow tumbling parameters.

This work could supplement further investigations into taming of active
turbulence. Moreover, it provides numerical context to both the work in Maitra
_et al._ (2018) as well as possible experimental assay results in which the
quadrupolar force plays a role in the active stabilization or destabilization
of the system. The exotic defect organisation phases revealed by the current
numerical study can further trigger new experimental exploration of the
topological defect structures in strongly confined active nematic systems.

## Acknowledgements

M. R. N. acknowledges the support of the Clarendon Fund Scholarship. A. D.
acknowledges support from the Novo Nordisk Foundation (grant No.
NNF18SA0035142), Villum Fonden (Grant no. 29476), and funding from the
European Union’s Horizon 2020 research and innovation program under the Marie
Sklodowska-Curie grant agreement No. 847523 (INTERACTIONS).

## References

  * Sanchez  _et al._ (2012) T. Sanchez, D. T. Chen, S. J. DeCamp, M. Heymann, and Z. Dogic, Nature 491, 431 (2012).
  * Weirich  _et al._ (2019) K. L. Weirich, K. Dasbiswas, T. A. Witten, S. Vaikuntanathan, and M. L. Gardel, Proceedings of the National Academy of Sciences 116, 11125 (2019).
  * Marchetti  _et al._ (2013) M. C. Marchetti, J. F. Joanny, S. Ramaswamy, T. B. Liverpool, J. Prost, M. Rao, and R. A. Simha, Reviews of Modern Physics 85, 1143 (2013).
  * Needleman and Dogic (2017) D. Needleman and Z. Dogic, Nature Reviews Materials 2, 1 (2017).
  * Jülicher  _et al._ (2018) F. Jülicher, S. W. Grill, and G. Salbreux, Reports on Progress in Physics 81, 076601 (2018).
  * Toner and Tu (1998) J. Toner and Y. Tu, Physical Review E 58, 4828 (1998).
  * Bate  _et al._ (2019) T. E. Bate, E. J. Jarvis, M. E. Varney, and K.-T. Wu, Soft Matter 15, 5006 (2019).
  * Duclos  _et al._ (2018) G. Duclos, C. Blanch-Mercader, V. Yashunsky, G. Salbreux, J.-F. Joanny, J. Prost, and P. Silberzan, Nature physics 14, 728 (2018).
  * Kawaguchi  _et al._ (2017) K. Kawaguchi, R. Kageyama, and M. Sano, Nature 545, 327 (2017).
  * Poujade  _et al._ (2007) M. Poujade, E. Grasland-Mongrain, A. Hertzog, J. Jouanneau, P. Chavrier, B. Ladoux, A. Buguin, and P. Silberzan, Proceedings of the National Academy of Sciences 104, 15988 (2007).
  * Meacock  _et al._ (2020) O. J. Meacock, A. Doostmohammadi, K. R. Foster, J. M. Yeomans, and W. M. Durham, Nature Physics 17, 205 (2020).
  * Wu  _et al._ (2017) K.-T. Wu, J. B. Hishamunda, D. T. N. Chen, S. J. DeCamp, Y.-W. Chang, A. Fernández-Nieves, S. Fraden, and Z. Dogic, Science 355 (2017), 10.1126/science.aal1979.
  * Aditi Simha and Ramaswamy (2002) R. Aditi Simha and S. Ramaswamy, Physical Review Letters 89, 058101 (2002).
  * Voituriez  _et al._ (2005) R. Voituriez, J.-F. Joanny, and J. Prost, EPL (Europhysics Letters) 70, 404 (2005).
  * Thampi  _et al._ (2015) S. P. Thampi, A. Doostmohammadi, R. Golestanian, and J. M. Yeomans, EPL (Europhysics Letters) 112, 28004 (2015).
  * Giordano  _et al._ (2021) M. G. Giordano, F. Bonelli, L. N. Carenza, G. Gonnella, and G. Negro, EPL (Europhysics Letters) 133, 58004 (2021).
  * Giomi  _et al._ (2011) L. Giomi, L. Mahadevan, B. Chakraborty, and M. F. Hagan, Physical Review Letters 106, 218101 (2011).
  * Thampi and Yeomans (2016) S. P. Thampi and J. M. Yeomans, The European Physical Journal Special Topics 225, 651 (2016).
  * Doostmohammadi  _et al._ (2018) A. Doostmohammadi, J. Ignés-Mullol, J. M. Yeomans, and F. Sagués, Nature communications 9, 1 (2018).
  * Brusatin  _et al._ (2018) G. Brusatin, T. Panciera, A. Gandin, A. Citron, and S. Piccolo, Nature materials 17, 1063 (2018).
  * Doostmohammadi and Yeomans (2019) A. Doostmohammadi and J. M. Yeomans, The European Physical Journal Special Topics 227, 2401 (2019).
  * Shankar  _et al._ (2021) S. Shankar, A. Souslov, M. J. Bowick, M. C. Marchetti, and V. Vitelli, “Topological active matter,” (2021), arXiv:2010.00364 [cond-mat.soft] .
  * Keber  _et al._ (2014) F. C. Keber, E. Loiseau, T. Sanchez, S. J. DeCamp, L. Giomi, M. J. Bowick, M. C. Marchetti, Z. Dogic, and A. R. Bausch, Science 345, 1135 (2014).
  * Metselaar  _et al._ (2019) L. Metselaar, J. M. Yeomans, and A. Doostmohammadi, Physical Review Letters 123, 208001 (2019).
  * Carenza  _et al._ (2020) L. N. Carenza, G. Gonnella, D. Marenduzzo, and G. Negro, Physica A: Statistical Mechanics and its Applications 559, 125025 (2020).
  * Hardoüin  _et al._ (2019) J. Hardoüin, R. Hughes, A. Doostmohammadi, J. Laurent, T. Lopez-Leon, J. M. Yeomans, J. Ignés-Mullol, and F. Sagués, Communications Physics 2, 1 (2019).
  * Wioland  _et al._ (2016) H. Wioland, E. Lushi, and R. E. Goldstein, New Journal of Physics 18, 075002 (2016).
  * Doostmohammadi  _et al._ (2016) A. Doostmohammadi, M. F. Adamer, S. P. Thampi, and J. M. Yeomans, Nat. Comm. 7, 10557 (2016).
  * Thijssen  _et al._ (2020) K. Thijssen, M. R. Nejad, and J. M. Yeomans, Physical Review Letters 125, 218004 (2020).
  * Mak  _et al._ (2016) M. Mak, M. H. Zaman, R. D. Kamm, and T. Kim, Nature Communications 7 (2016).
  * Ramaswamy (2010) S. Ramaswamy, Annual Review Condensed Matter Physics 1, 323 (2010).
  * Prost  _et al._ (2015) J. Prost, F. Jülicher, and J.-F. Joanny, Nature Physics 11, 111 (2015).
  * Mathijssen  _et al._ (2016) A. J. T. M. Mathijssen, A. Doostmohammadi, J. M. Yeomans, and T. N. Shendruk, Journal of Fluid Mechanics 806, 35–70 (2016).
  * Maitra  _et al._ (2018) A. Maitra, P. Srivastava, M. C. Marchetti, J. S. Lintuvuori, S. Ramaswamy, and M. Lenz, Proceedings of the National Academy of Sciences 115, 6934 (2018).
  * Saw  _et al._ (2017) T. B. Saw, A. Doostmohammadi, V. Nier, L. Kocgozlu, S. Thampi, Y. Toyama, P. Marcq, C. T. Lim, J. M. Yeomans, and B. Ladoux, Nature 544, 212 (2017).
  * Copenhagen  _et al._ (2020) K. Copenhagen, R. Alert, N. S. Wingreen, and J. W. Shaevitz, Nature Physics 17, 211 (2020).
  * Maroudas-Sacks  _et al._ (2020) Y. Maroudas-Sacks, L. Garion, L. Shani-Zerbib, A. Livshits, E. Braun, and K. Keren, Nature Physics 17, 251 (2020).
  * Zhang  _et al._ (2021) J. Zhang, N. Yang, P. K. Kreeger, and J. Notbohm, APL Bioengineering 5, 036103 (2021).
  * Blanch-Mercader  _et al._ (2021) C. Blanch-Mercader, P. Guillamat, A. Roux, and K. Kruse, Physical Review Letters 126, 028101 (2021).
  * Doostmohammadi and Ladoux (2021) A. Doostmohammadi and B. Ladoux, Trends in Cell Biology (2021).
  * Koch and Subramanian (2011) D. L. Koch and G. Subramanian, Annual Review of Fluid Mechanics 43, 637 (2011).
  * Kumar  _et al._ (2018) N. Kumar, R. Zhang, J. J. de Pablo, and M. L. Gardel, Science Advances 4 (2018).
  * Turiv  _et al._ (2020) T. Turiv, R. Koizumi, K. Thijssen, M. M. Genkin, H. Yu, C. Peng, Q.-H. Wei, J. M. Yeomans, I. S. Aranson, A. Doostmohammadi, and O. D. Lavrentovich, Nature Physics 16, 481 (2020).
  * Bäuerle  _et al._ (2020) T. Bäuerle, R. C. Löffler, and C. Bechinger, Nature Communications 11 (2020).
  * Scholz  _et al._ (2018) C. Scholz, M. Engel, and T. Pöschel, Nature Communications 9, 1497 (2018).
  * Chaté (2020) H. Chaté, Annual Review of Condensed Matter Physics 11, 189 (2020).
  * Junot  _et al._ (2017) G. Junot, G. Briand, R. Ledesma-Alonso, and O. Dauchot, Physical Review Letters 119, 028002 (2017).
  * Kumar  _et al._ (2014) N. Kumar, H. Soni, S. Ramaswamy, and A. Sood, Nature Communications 5, 4688 (2014).
  * Peruani  _et al._ (2012) F. Peruani, J. Starruß, V. Jakovljevic, L. Søgaard-Andersen, A. Deutsch, and M. Bär, Physical Review Letters 108, 098102 (2012).
  * Großmann  _et al._ (2016) R. Großmann, F. Peruani, and M. Bär, Physical Review E 94, 050602 (2016).
  * Nagai (2018) K. H. Nagai, Biophysics and Physicobiology 15, 51 (2018).
  * Alvarado  _et al._ (2014) J. Alvarado, B. M. Mulder, and G. H. Koenderink, Soft Matter 10, 2354 (2014).
  * Duclos  _et al._ (2017) G. Duclos, C. Erlenkämper, J.-F. Joanny, and P. Silberzan, Nature Physics 13, 58 (2017).
  * Mueller  _et al._ (2019) R. Mueller, J. M. Yeomans, and A. Doostmohammadi, Physical Review Letters 122, 048004 (2019).
  * Blanch-Mercader  _et al._ (2018) C. Blanch-Mercader, V. Yashunsky, S. Garcia, G. Duclos, L. Giomi, and P. Silberzan, Physical Review Letters 120, 208101 (2018).
  * Oza and Dunkel (2016) A. U. Oza and J. Dunkel, New Journal of Physics 18, 093006 (2016).
  * Putzig  _et al._ (2016) E. Putzig, G. S. Redner, A. Baskaran, and A. Baskaran, Soft matter 12, 3854 (2016).
  * Ternet  _et al._ (1999) D. J. Ternet, R. G. Larson, and L. G. Leal, Rheologica Acta 38, 183 (1999).
  * Beris and Edwards (1994) A. N. Beris and B. J. Edwards,  _Thermodynamics of Flowing Systems_ (Oxford University Press, 1994).
  * Nejad  _et al._ (2021) M. R. Nejad, A. Doostmohammadi, and J. M. Yeomans, Soft Matter 17, 2500 (2021).
  * Thampi  _et al._ (2014) S. P. Thampi, R. Golestanian, and J. M. Yeomans, EPL (Europhysics Letters) 105, 18001 (2014).
  * (62) “Supplementary Movies,” https://www.dropbox.com/sh/alwauc0joang94r/AAC2_Ue1Db9W2e1wx406Iwzga?dl=0.
  * Blow  _et al._ (2014) M. L. Blow, S. P. Thampi, and J. M. Yeomans, Physical Review Letters 113, 248303 (2014).
  * Pearce  _et al._ (2021) D. Pearce, J. Nambisan, P. Ellis, A. Fernandez-Nieves, and L. Giomi, Physical Review Letters 127, 197801 (2021).
  * De Gennes and Prost (1993) P.-G. De Gennes and J. Prost,  _The physics of liquid crystals_ , Vol. 83 (Oxford university press, 1993).
  * Shankar and Marchetti (2019) S. Shankar and M. C. Marchetti, Physical Review X 9, 041047 (2019).
  * Li  _et al._ (2021) Z.-Y. Li, D.-Q. Zhang, and B. Li, Physical Review Research 3, 023253 (2021).
  * Shankar  _et al._ (2018) S. Shankar, S. Ramaswamy, M. C. Marchetti, and M. J. Bowick, Physical Review Letters 121, 108002 (2018).