# Active smectics on a sphere

Michael Nestler1, Simon Praetorius1, Zhi-Feng Huang2, Hartmut Löwen3, Axel
Voigt1,4 1 Institute of Scientific Computing, Technische Universität Dresden,
01062 Dresden, Germany
2 Department of Physics and Astronomy, Wayne State University, Detroit,
Michigan 48201, USA
3 Institut für Theoretische Physik II: Weiche Materie, Heinrich-Heine-
Universität Düsseldorf, 40225 Düsseldorf, Germany
4 Center for Systems Biology Dresden, Pfotenhauerstr. 108, 01307 Dresden,
Germany<EMAIL_ADDRESS>
###### Abstract

The dynamics of active smectic liquid crystals confined on a spherical surface
is explored through an active phase field crystal model. Starting from an
initially randomly perturbed isotropic phase, several types of topological
defects are spontaneously formed, and then annihilate during a coarsening
process until a steady state is achieved. The coarsening process is highly
complex involving several scaling laws of defect densities as a function of
time where different dynamical exponents can be identified. In general the
exponent for the final stage towards the steady state is significantly larger
than that in the passive and in the planar case, i.e., the coarsening is
getting accelerated both by activity and by the topological and geometrical
properties of the sphere. A defect type characteristic for this active system
is a rotating spiral of evolving smectic layering lines. On a sphere this
defect type also determines the steady state. Our results can in principle be
confirmed by dense systems of synthetic or biological active particles.

##  1 Introduction

Liquid crystals consist of particles that can possess both orientational and
positional degrees of freedom [1, 2]. This induces a wealth of “mesophases”
with partial orientational and positional ordering. In general these phases
contain topological defects. For pure positional ordering such defects are
characterized by deviations from the equilibrium number of neighbors, and for
orientational order alone a director field can be defined and topological
defects emerge where the orientation of this field shows a singularity. In
smectic phases, orientational and positional orders are present simultaneously
and compete with each other, and these defects also interact. Describing and
predicting the coarsening of defected, multidomain smectic configurations is a
long-lasting problem. Similar to grain growth in crystalline materials, self-
similarity in domain coarsening has been found, which provides a scaling law
$t^{\alpha}$ for the characteristic domain size with time $t$ and a scaling
exponent $\alpha$. The exponent $\alpha$ varies between $1/5$ and $1/2$
according to different theoretical [3, 4, 5, 6, 7, 8] and experimental [9, 10]
investigations in 2D space. The coarsening process is usually dominated by the
motion of grain boundaries, which move over large distances and are driven by
the curvature of the smectic lines and internal distortions within domains.
These distortions are positional defects, in particular disclinations and
dislocations. The motion of dislocations [11] and the motion of disclinations
[12, 8] in smectics are coupled by topological constraints and their
interactions influence the coarsening process. This complexity might explain
the broad spectrum of reported values of scaling exponent. Results on power
law scaling for defect densities in 2D space are more narrow, for which the
disclination density scales as $t^{-1/2}$ and the dislocation density scales
as $t^{-1/3}$ (see [9, 12, 7, 8]).

Here we address two questions: How is this coarsening process influenced if
the smectic phase forms on the surface of a sphere, with the director
constrained to be tangential, and how does it change if the smectic phase is
active? The first introduces topological constraints resulting from the
Poincaré-Hopf theorem and Euler’s theorem on polyhedra, and forces the defects
(orientational and positional) to be present also in equilibrium. The second
brings the system into nonequilibrium and there is no relaxation towards
equilibrium but rather towards a nonequilibrium steady state. Smectics on
spherical surfaces, the so-called smectic shells, have been experimentally
realized [13, 14]. The topological and geometrical frustration in these
experiments indeed leads to the formation of defects. However, in these
experiments the thickness of the smectic layer and the proposed anchoring
conditions for the director also influence the equilibrium state. Our study
here corresponds to the vanishing thickness limit of an ideal smectic sphere.
This case is characterized by two $+1$ (orientational) defects located on
opposing poles of the sphere and a director field aligned with geodesic lines
connecting these defects. Activity in smectics can be realized by Janus rods
which are, for instance, catalytically driven and thus convert chemical energy
into mechanical work [15]. Realizations of active smectics are mostly
considered in 2D flat space [16, 17]. These systems can be characterized by
aligning interaction [18], nonreciprocality [19], or nonlinear mutual feedback
[20]. Due to the nontrivial coupling between orientational and positional
degrees of freedom with the active driving force, active smectic systems and
their defect dynamics are more complex than other active liquid crystals,
e.g., active nematics [21]. To the best of our knowledge, studies of defects
in active smectics if constrained to the surface of a sphere, are still
lacking. Potential applications for active smectics on surfaces are found in
biology. As mentioned in [22], developing tissues on substrates or myosin
filaments in cellular cortex are important examples. The investigation of such
systems on a sphere is the first step to understand the complex interplay of
orientational and positional defects in multidomain smectic systems under
topological constraints, geometric effects, and active driving.

We study this complex interplay computationally using an active phase field
crystal (PFC) model on a sphere [23]. In the following we first review the
known results for equilibrium configurations on a sphere and summarize results
for the dynamics of the corresponding active systems. We then introduce the
active PFC model which has recently been used to model active smectics in 2D
space [24], reformulate it on surfaces and simulate the coarsening process on
a sphere. While essential processes during coarsening are explained on highly
resolved simulations, results on scaling behavior are addressed by large scale
statistical investigations. To unveil the influence of the spherical topology
and activity, we compare the results with corresponding simulations in flat 2D
space and vary the strength of activity.

##  2 Review of known results

###  2.1 Equilibrium configurations on a sphere

Equilibrium configurations for pure orientational or pure positional ordering
on a sphere have been intensely studied (see, e.g., [25] for a review). To
identify the optimal distribution of interacting particles on a sphere is a
classical problem [26, 27, 28], for which above a certain number of particles
there exist non-trivial configurations of defects, the so-called grain
boundary scars [29, 30, 31], all consistent with Euler’s theorem on
polyhedrons. For anisotropic particles the situation becomes more complicated.
For orientational ordering [32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42],
typically four $+1/2$ defects arrange in a tetrahedral configuration,
maximizing their distance from each other. Much less is known for structures
that are both orientationally and positionally ordered, such as the smectic
phase in curved space. Clearly the simultaneous presence and competition of
both local orientational and positional orders put this situation into a much
more complex category. For results on a sphere we refer to [36, 43, 44, 45,
46]. Particle-resolved computer simulations of hard rods indicate a wealth of
different states, including equatorial smectics with two $+1$ defects at the
poles, consistent with the finding in [13, 14], and the configurations with
four $+1/2$ or six $+1/2$ and two $-1/2$ defects. In these situations the
strain at the poles is relieved by separating closely positioned half-strength
defects. All these configurations obviously fulfill the topological constraint
of the Poincaré-Hopf theorem.

###  2.2 Active liquid crystals on a sphere

The active liquid crystalline phases, which have been of tremendous interest
in recent years [47, 15, 48, 49, 50], exhibit qualitatively different
behaviors as compared to those in the corresponding passive systems. Their
governing dynamics is non-relaxational, as in the far-from-equilibrium pattern
formation processes [51]. This type of system on curved surfaces has been
termed topological active matter [52] as it combines active liquid crystals
with topological constraints. Again, studies are mainly concerned with either
pure orientational order, i.e., active nematics [52, 53, 54, 55, 56] showing
persistent oscillations between the tetrahedral and a planar defect
configuration as one possible dynamical process, or the appearance of polar
order which leads to polar vortex and circulating band states [57, 58], or
pure positional order as in active crystals [23], which can lead to collective
rotation as seen in epithelia tissue confined on a sphere [59, 60, 61, 62].
While in flat space activity can lead to continuous creation and annihiliation
of defects [21], in topological active matter such creation of defects is
suppressed over a wide range of activity.

###  2.3 Phase field crystal and other modeling approaches

In [63] a phase field crystal (PFC) model for liquid crystals was introduced,
and due to various couplings between the rescaled density field, the local
nematic order parameter, and the mean local direction of the orientations and
their gradients, various stable phases of the equilibrium free-energy
functional could be identified, including the smectic-A phase. Formulating
this model on a sphere allows us to analyse these phases under the topological
constraint and the influence of curvature. Besides the equilibrium
configurations, the dynamics of these liquid crystalline phases is of
importance, which could be addressed by the diffusive dynamics of the PFC
model for liquid crystals. However, results mainly exist for systems with pure
orientational or positional order, ranging from the classical Lifshitz-
Slyozov-Wagner theory [64], to large scale PFC simulations [65, 66] and
simulations for nematic liquid crystals on surfaces [41, 42]. They are all
characterized by annihilation of defects and convergence towards a (local)
minimum with energetically prefered defects and the corresponding spatial
arrangement. In addition, active PFC models [67, 68, 69, 70, 71, 72, 73, 24]
have been introduced. These models are related to the PFC model for liquid
crystals [63] and have been formulated on a sphere in [74, 23] to study active
crystals.

##  3 Results

In the following we first introduce the surface active PFC model and then
describe large-scale simulation results for active smectics on a sphere. We
compare results of defect dynamics in multidomain smectic configurations,
their coarsening and the emerging of stable configurations to the situation in
flat space with periodic boundary conditions.

###  3.1 Model

In [23] the active PFC model introduced in [67] has been formulated on a
sphere ${\mathcal{S}}=R\,{\mathbb{S}_{2}}$ with radius $R$, where
${\mathbb{S}_{2}}$ is the two dimensional unit sphere embedded in
$\mathbb{R}^{3}$. The model was used to study active crystals. Here, we
consider a different parameter setting which is applicable for active
smectics. In the model the position $\mathbf{r}$ on the sphere ${\mathcal{S}}$
is parametrized by $\mathbf{r}(\theta,\phi)=R\mathbf{\hat{u}}(\theta,\phi)$
with the orientational unit vector
$\mathbf{\hat{u}}(\theta,\phi)=(\sin\theta\cos\phi,\sin\theta\sin\phi,\cos\theta)^{\mathrm{T}}$
and the spherical coordinates $\theta\in[0,\pi]$ and $\phi\in[0,2\pi)$. We
consider a rescaled density field $\psi(\mathbf{r},t)$ and a polarization
field $\mathbf{p}(\mathbf{r},t)$, which is tangential to ${\mathcal{S}}$ at
$\mathbf{r}$, i.e.,
$\mathbf{p}(\mathbf{r},t)=p_{\theta}(\mathbf{r},t)\partial_{\theta}\mathbf{\hat{u}}+p_{\phi}(\mathbf{r},t)\partial_{\phi}\mathbf{\hat{u}}\in{{\mathrm{T}_{\mathbf{r}}{{\mathcal{S}}}}}$
with scalar polarization fields $p_{\theta}(\mathbf{r},t)$ and
$p_{\phi}(\mathbf{r},t)$ and the tangent space
${{\mathrm{T}_{\mathbf{r}}{{\mathcal{S}}}}}$ of the sphere ${\mathcal{S}}$ at
point $\mathbf{r}$.

The free-energy functionals read

$\displaystyle\mathcal{F}_{\psi}[\psi]$
$\displaystyle=\int_{{\mathcal{S}}}\\!\\!\Big{\\{}\frac{1}{2}\psi\big{[}{\varepsilon}+(q_{0}+\boldsymbol{\triangle}_{{\mathcal{S}}})^{2}\big{]}\psi+\frac{1}{4}\psi^{4}\Big{\\}}\mathrm{d}\mathcal{S}\,,$
(1) $\displaystyle\mathcal{F}_{\mathbf{p}}[\mathbf{p}]$
$\displaystyle=\int_{{\mathcal{S}}}\\!\\!\Big{(}\frac{1}{2}C_{1}\lVert\mathbf{p}\rVert^{2}+\frac{1}{4}C_{2}\lVert\mathbf{p}\rVert^{4}\Big{)}\mathrm{d}\mathcal{S}\,,$
(2)

with $\epsilon<0$, the characteristic wave number $q_{0}=1$ after rescaling,
$C_{1}>0$ which leads to the suppression of any spontaneous ordering of
orientational alignment, and $C_{2}$ modeling the higher order nonlinear term
which is neglected in the following. The gradient and divergence operators are
considered in spherical coordinates, i.e.,

$\displaystyle{\operatorname{grad}_{{\mathcal{S}}}}\psi$
$\displaystyle=\frac{1}{R}\Big{[}(\partial_{\theta}\mathbf{\hat{u}})\partial_{\theta}\psi+\frac{1}{\sin^{2}\theta}(\partial_{\phi}\mathbf{\hat{u}})\partial_{\phi}\psi\Big{]},$
(3) $\displaystyle{\operatorname{div}_{{\mathcal{S}}}}\mathbf{p}$
$\displaystyle=\frac{1}{R}\big{(}\cot\theta
p_{\theta}+\partial_{\theta}p_{\theta}+\partial_{\phi}p_{\phi}\big{)},$ (4)

and the surface Laplace-Beltrami operator is defined as
$\boldsymbol{\triangle}_{{\mathcal{S}}}={\operatorname{div}_{{\mathcal{S}}}}\,{\operatorname{grad}_{{\mathcal{S}}}}$.

The dynamic equations describe the active-particle transport tangential to
${\mathcal{S}}$ and read

$\displaystyle\partial_{t}\psi$
$\displaystyle=\boldsymbol{\triangle}_{{\mathcal{S}}}\frac{\delta\mathcal{F}_{\psi}}{\delta\psi}-v_{0}\,{\operatorname{div}_{{\mathcal{S}}}}\mathbf{p}\,,$
(5) $\displaystyle\partial_{t}\mathbf{p}$
$\displaystyle=-(\boldsymbol{\triangle}_{\mathrm{dR}}+D_{r})\frac{\delta\mathcal{F}_{\mathbf{p}}}{\delta\mathbf{p}}-v_{0}\,{\operatorname{grad}_{{\mathcal{S}}}}\psi\,,$
(6)

with $v_{0}$ the strength of particle self-propulsion and $D_{r}$ the
rotational diffusion constant. The vector-Laplacian
$\boldsymbol{\triangle}_{\mathrm{dR}}$ is the surface Laplace-de Rham operator
$\boldsymbol{\triangle}_{\mathrm{dR}}=-{\operatorname{grad}_{{\mathcal{S}}}}{\operatorname{div}_{{\mathcal{S}}}}-{\operatorname{rot}_{{\mathcal{S}}}}{\operatorname{Rot}_{{\mathcal{S}}}}$,
where

$\displaystyle{\operatorname{rot}_{{\mathcal{S}}}}\psi$
$\displaystyle=\frac{1}{R\sin\theta}\big{[}-(\partial_{\theta}\mathbf{\hat{u}})\partial_{\phi}\psi+(\partial_{\phi}\mathbf{\hat{u}})\partial_{\theta}\psi\big{]},$
(7) $\displaystyle{\operatorname{Rot}_{{\mathcal{S}}}}\mathbf{p}$
$\displaystyle=\frac{1}{R}\Big{(}2\cos\theta
p_{\phi}\\!-\\!\frac{1}{\sin\theta}\partial_{\phi}p_{\theta}+\sin\theta\partial_{\theta}p_{\phi}\Big{)},$
(8)

are the surface curl operators in spherical coordinates. The equations have
been rescaled with a diffusive timescale and a length scale set via the
pattern periodicity, and in this form constitute the minimal continuum field
model which was used in [23] to model active crystals on a sphere. When
$v_{0}=0$, Eq. (5) reduces to the surface PFC model, which has been applied to
describe optimal ordering and crystallization of passive particles on a sphere
[31, 75]. Equations (5) and (6) reduces to the active PFC model [67] in flat
space in the limit $R\to\infty$. Within this limit, the model was used in [24]
to simulate defect dynamics in two-dimensional active smectics. This has been
done by considering a parameter setting of $\epsilon=-0.98$ and the average
density $\psi_{0}=0$ which favours a travelling stripe phase [67]. Such a
stripe phase can also emerge on a sphere. Regarding the remaining model
parameters we follow the aforementioned studies [76, 23] and choose
$D_{r}=0.5$, $C_{1}=0.2$, and $C_{2}=0$. For spherical geometry we consider a
sphere with radius $R=100$, such that for passive case $50$ repetitions of a
rotational symmetric smectic pattern perfectly fit to the surface. In the flat
case a square periodic domain is simulated with size $L$, where $L$ is chosen
to match domain area w.r.t the sphere.

We use vector spherical harmonics, details in the methods section 5, to
conduct multiple simulations ($15$ for flat geometry and $50$ for spherical
geometry) of the coarsening process for each considered activity
$v_{0}\in\\{0.305,\,0.35,\,0.425,\,0.5\\}$. Each simulation starts from
uniform $\psi(0)=\psi_{0}$ disturbed by random noise $\psi^{\prime}$ with
magnitude $|\psi^{\prime}|<0.01$ applied at each grid point.

In order to identify orientational defects we define a nematic order parameter
almost everywhere $\mathbf{x}\mbox{ :
}\|{\operatorname{grad}_{{\mathcal{S}}}}\psi(\mathbf{x})\|\neq 0$:

$\displaystyle\mathbf{Q}(\mathbf{x})=\mathbf{d}(\mathbf{x})\mathbf{d}(\mathbf{x})-\frac{1}{2}\boldsymbol{\Pi}(\mathbf{x}),\quad\mathbf{d}(\mathbf{x})=\frac{{\operatorname{grad}_{{\mathcal{S}}}}\psi(\mathbf{x})\times\boldsymbol{\nu}(\mathbf{x})}{\|{\operatorname{grad}_{{\mathcal{S}}}}\psi(\mathbf{x})\|},$
(9)

with
$\boldsymbol{\Pi}(\mathbf{x})=\mathbf{I}-\boldsymbol{\nu}(\mathbf{x})\boldsymbol{\nu}(\mathbf{x})$
the surface identity of the considered geometry and $\boldsymbol{\nu}$ the
outward oriented surface normal. The nematic texture is then extended to the
complete domain by smooth continuation. An associated distortion energy can be
defined along the tensor components $I,J$ of $\mathbf{Q}$ w.r.t the
coordinates of embedding space $\mathbb{R}^{3}$

$\displaystyle{\mathcal{F}}_{\mathbf{Q}}[\mathbf{Q}]=\int_{{\mathcal{S}}}f_{d}(\mathbf{Q})\mathrm{d}\mathcal{S},\quad
f_{d}(\mathbf{Q})=\frac{1}{2}\sum_{I,J=1}^{3}\|{\operatorname{grad}_{{\mathcal{S}}}}\mathbf{Q}_{I,J}\|^{2},$
(10)

which will be used to identify defects. For details we refer to the methods
section 5.

###  3.2 Pattern formation in active smectics

In order to discuss the complex dynamics of smectic pattern formation and
coarsening we introduce the following notions:

  * •

Smectic lines are lines of minimal (maximal) density $\psi$. In equilibrium
they are quasi parallel and equidistant. In the active system these lines move
consistently in a direction orthogonal to themselves;

  * •

Grains are regions with quasi parallel and equidistant smectic lines. They are
characterized by the orientation of the lines or the direction of their
movement;

  * •

Grain boundaries separate grains with different line orientation or movement
direction from each other. If more than two grains are in contact the grain
boundaries meet in junction points;

  * •

Defects are localized areas where the smectic lines are not parallel or not
equidistant. For each defect a topological charge is associated;

  * •

Defect clusters are compact aggregates of defects. In contrast to grain
boundaries, which are typically chainlike aggregations of defects, these
clusters lead to rotation;

  * •

Rotating grains are a combination of defect clusters and their associated
grains. Due to the rotation of the cluster, smectic lines form an Archimedes
spiral with two rotation directions and two movement directions of the smectic
lines. They move either inward or outward.

Fig. 1 and 2 highlight most of these features. All phenomena introduced above
can be observed in the video corresponding to Fig. 2 (see SI).

Figure 1: Smectic pattern formation and dynamics in sphere domain. [A]
Snapshots of $\psi$ exhibiting a typical set of phenomena with [A1] a single
disclination in smectic line structure with topological charge of $+1/2$
(location denoted by black semi circle); [A2] a grain boundary, where the
nematic texture exhibits an alternating pattern of $+1/2$ (black semi circles)
and $-1/2$ (black tribars) defects; [A3] an isolated dislocation in smectic
lines, for which the associated nematic texture consists of two disclination
defects with opposite sign yielding an effective topological charge $0$ for
the dislocation. [B] A rotating cluster formed by several disclinations, where
panels 1–4 depict the temporal evolution of $\psi$ with time steps of
$\triangle t=100$, showing an outward spiraling movement of lines (marked in
green and pink) with rearranging disclinations and the associated nematic
defects at center. The total topological charge of the cluster remains as $+1$
and its center does not move.

Fig. 2 and the corresponding video in SI show the coarsening process for one
simulation with $v_{0}=0.5$. Within Fig. 2 at $t_{0}=500$ smectic lines have
emerged and formed multiple grains. These grains move orthogonal to the
smectic lines, with black arrows indicating the direction of movement. These
phenomena are similar to the corresponding state in flat 2D space and are also
reminiscent to the known coarsening processes in passive systems. At
$t_{1}=900$, at the front side (top panel) grain boundaries curl up and form
rotating grains which induce outward movement w.r.t the central rotating
defect cluster. On the back side (bottom panel) the direction of the grain
movement changes without the curling up of grain boundaries. At $t_{2}=1250$
two counterclockwise rotating grains have formed with outward moving smectic
lines, inducing a long grain boundary between these two incompatible grains.
Note that the rotating defect clusters do not move as soon as they are formed.
After a metastable state where grain boundaries wiggle, at $t_{3}=4000$ the
rotating grain with its center at the back side expands and destroys the
rotating cluster at the front side. Snapshot at $t_{4}=9500$ depicts an
intermediate configuration where the front cluster is destroyed, yet the
remaining grain moves in the clockwise direction and collapses. This process
leads to the final single grain configuration at $t_{5}=20000$ consisting of a
counterclockwise rotating cluster at the back side with outward moving smectic
lines. On the opposite side, with maximum geodesic distance, a clockwise
rotating cluster has formed with inward moving smectic lines. The
corresponding transient of the distortion energy ${\mathcal{F}}_{\mathbf{Q}}$,
together with the considered time instances, is shown in Fig. 3 (top panel of
Simulation 4).

Figure 2: Grain growth during the coarsening process for $v_{0}=0.5$ in sphere
domain, showing as a sequence of snapshots depicting $\psi$ at various times
marked in Simulation 4 of Fig. 3, with front view (top panels) and back view
(bottom) in terms of $180^{\circ}$ rotation w.r.t vertical axis. Black arrows
indicate the direction of grain movement and black circle arrows mark the
rotating defect clusters. Black bold letters (A, B, C) indicate rotating
clusters and correspond to the marks given in Fig. 4 [B]. Snapshots are taken
at $t_{0}=500$, $t_{1}=900$, $t_{2}=1250$, $t_{3}=4000$, $t_{4}=9500$, and
$t_{5}=20000$. An animation of the presented snapshots and intermediate grain
configurations is provided in the SI.

###  3.3 Distinct coarsening regimes

Figure 3: Sample coarsening processes for $v_{0}=0.5$: In sphere (top) and
plane (bottom) domain, showing the time evolution of nematic distortion energy
$\mathcal{F}_{\mathbf{Q}}(t)$ for four different simulations. Diamond markers
in Simulation 4 correspond to the snapshots used in Fig. 2 and Fig. 4. Labels
a, b, c and d mark distinct coarsening regimes. The scaling exponent for the
plane domain is roughly $-1/2$ in all four simulations, as indicated by the
dashed lines. Figure 4: Dynamics of defect distribution during the coarsening
process: Left panel [A] for planar and right panel [B] for sphere domain for
$v_{0}=0.5$. The different time instances $t_{0}$, $\cdots$, $t_{5}$
correspond to the transients marked in Fig. 3 (Simulation 4) showing the
nematic distortion energy $\mathcal{F}_{\mathbf{Q}}$. Dark blue patches
indicate smectic lines that are free from nematic defects, while distinct
grain boundaries are indicated by green to yellow lines. Movements of isolated
defects are visible as short bright strokes. Rotating defect clusters and
directions of defect motion are marked by white arrows. A video of the
corresponding time evolution is provided in the SI.

The transients of the distortion energy ${\mathcal{F}}_{\mathbf{Q}}$ in the
spherical case, as shown in Fig. 3 (top panels), indicate different coarsening
regimes. While different simulation runs differ quantitatively, they share
common qualitative features. We can identify four distinct regimes which are
marked in each plot by labels a, b, c and d. Regime a refers to initial
coarsening, as characterized by small grains which merge or dissolve. Defect
clusters are rarely observed in this regime. The scaling exponent in this
regime is roughly $-1/2$ in all simulations. In regime b, a set of large
rotating clusters is formed, while small grains without rotating clusters are
rarely present. The spatial arrangement of the clusters and their rotation
direction seem random but impact the subsequent evolution. This regime is
characterized by a coarsening of rotating clusters. As in Ostwald ripening,
large rotating clusters grow on the expense of small rotating clusters. This
leads to a faster coarsening compared to regime a. The associated rotating
grains are more stable than grains without rotating clusters. A critical state
is reached when a grain boundary reaches a rotating cluster and dissolves it.
The remaining grain is then typically dissolved into the neighbouring rotating
grain. Some intermediate configurations during the transition from regime a to
regime b are more stable or long-lasting than others. In the configuration
shown in Fig. 2, after a similar coarsening process indicated by b1 in Fig. 3,
two almost equally sized rotating grains cover the domain, although their
orientations are different, with a long grain boundary separating them. Due to
their similar size the grain boundary migrates slowly, yielding a long time
period with quasi stable configuration, as indicated by b2 in Fig. 3. In
regime c, after the motion of grain boundary of the largest rotating grain has
dissolved any other rotating clusters, only a single grain and the
corresponding rotating cluster remain. The grain boundary then quickly
collapses, and an inward-directed rotating defect cluster appears on the
opposite pole position of the dominant rotating defect cluster. The collapse
of the grain boundary happens very fast in a zipping mode, leading to a large
coarsening exponent (see Section 3.4). The final regime d is featured by a
stable configuration of two rotating defect clusters (one inward, one outward)
positioned on opposing sides/poles of the sphere with maximized distance. The
resulting rotation axis is determined by the dominant rotating defect cluster
emerging during the transition from regime a to regime b. Also, the number of
defects in these final-state rotating defect clusters can vary in different
simulations and is determined by the formation of the dominant rotating defect
cluster.

In order to unveil the qualitative difference of the plane vs. sphere domain,
we conduct the same type of simulations on a planar 2D domain. The
corresponding distortion energy is shown in Fig. 3 (bottom panels). They do
not lead to distinct coarsening regimes and essentially can be described by a
single regime a, with scaling exponent roughly $-1/2$ in all simulations (as
indicated in the bottom panels), which terminates into the final configuration
d.

To further understand the differences we visualize the dynamics of defect
distributions in the temporally binned data around different time instances
for one sample simulation, as shown in Fig. 4, and compare the configurations
on the sphere and the plane domains. For the sphere the considered time
instances correspond to those of Fig. 2. These defect distributions are
obtained by collecting defect positions around the selected time instances
$t_{0},\ldots,t_{5}$, each of which corresponds to one of the $100$ temporal
bins (logarithmic sizes) set up between $t=10^{2}$ and $t=10^{5}$, and
evaluating a spatial distribution $D_{i}(x_{0},x_{1})$ ($100\times 100$
spatial bins) of the positions. The resulting two-dimensional probability
distribution of defect positions is plotted as log-scale colormap (dark blue
to yellow) with range limited to $[10^{-3},1]$. Dark blue patches indicate
smectic lines that are free from nematic defects, while distinct grain
boundaries are indicated by green to yellow lines. Movement of isolated
defects (disclinations or dislocations) is visible as short bright strokes.
Rotating clusters are marked and indicated by dashed circles. On the sphere
local coordinates $(x_{0},x_{1})$ are derived from spherical coordinates
($x_{0}$ corresponds to polar angle $\theta$, $x_{1}$ to azimuthal angle
$\phi$) such that strong distortions can appear near left and right boundaries
of snapshots.

Initially the situation is similar on the sphere and the plane domain. At the
early time $t_{0}$ an initial establishment of smectic lines in small patches
is visible for both domains. These patches merge and rearrange to larger
grains at time $t_{1}$. For the sphere domain, around time $t_{1}$ two
rotating clusters are formed, marked as A and B. This occurs by the curling of
grain boundaries. At a later time $t_{2}$ they evolve to two similarly sized
grains, yielding a metastable state where grain boundaries hardly move.
Eventually the grain boundaries move towards cluster A at $t_{3}$. In the
$t_{4}$ snapshot cluster A has been dissolved and the grain boundary collapses
with another rotating cluster. At $t_{5}$ the final minimal configuration
(ground state dictated by topology) of an outward rotating cluster B and an
inward rotating cluster C is formed. This evolution process corresponds to the
snapshots given in Fig. 2. Especially the collapse of the grain boundary at
$t_{4}$ is characteristic of the fast decay of $\mathcal{F}_{\mathbf{Q}}$
shown in Fig. 3.

The situation differs on the plane domain. At $t_{1}$ there are two grains
each containing a rotating cluster, labeled A and B. In transition from
$t_{1}$ to $t_{2}$ grains expand and merge, and only the grains containing
rotating clusters A and B remain. At $t_{3}$ the moving grain boundary reaches
a rotating cluster and dissolves it (as happened for cluster B in transition
from $t_{2}$ to $t_{3}$). With the vanishing of the last rotating cluster A
the grain boundaries dissolve into isolated defects and the collective motion
of smectic lines starts. In the $t_{4}$ snapshot parallel bright strokes
indicate collective movement of isolated defects along the direction indicated
by arrows. They rearrange and move collectively at the final configuration
shown at $t_{5}$. Identifying these details for the other simulation runs of
Fig. 3 leads to the same qualitative features. We never observe a collapsing
grain boundary in a planar 2D geometry as on the sphere domain at $t_{4}$. We
therefore associate this phenomenon of grain collapse with the topological and
geometrical constraints of the sphere.

###  3.4 Fast coarsening regime

We next examine how the fast coarsening regime identified above is influenced
by the activity and how it influences statistical coarsening laws. We
therefore conduct more simulations at different values of $v_{0}$, and again
compare our results with the corresponding simulations in flat space. Fig. 5
[A] shows the results of averaged distortion energy
${\mathcal{F}}_{\mathbf{Q}}$ over time, which can be related to the statistics
on the time decay of defect density.

Figure 5: Average behavior of coarsening for sphere and plane domains. [A] The
nematic distortion energy $\mathcal{F}_{\mathbf{Q}}$ averaged over ensemble of
simulation runs on sphere (blue lines, 50 runs for each activity) and plane
(red lines, 10 runs for each activity). Black dashed lines indicate the fitted
power-law decay. Black diamonds and dots indicate the fitting interval for the
fast coarsening regime and the transition to the steady state. [B] Time
evolution of $\mathcal{F}_{\mathbf{Q}}$ at $v_{0}=0.5$, for individual runs
(light color lines) and the corresponding average (solid black line). Left and
right panels correspond to sphere and plane geometry, respectively. [C]
Summary of values of scaling exponent $r$ obtained from the fitting in [A] as
a function of activity $v_{0}$, for sphere (blue) and plane (red) geometries.

For the planar geometry there is no clear distinction between different
scaling regimes, as shown in both Fig. 3 (bottom panels) and Fig. 5 [A]. This
is similar to passive smectics in 2D space. At small activity $v_{0}=0.305$ we
obtain the value of coarsening exponent of $0.26$, close to the scaling of
$t^{-1/3}$ for dislocations in 2D passive smectics [7, 8]. At large enough
activity the scaling is roughly $t^{-1/2}$.

For the sphere the initial coarsening regime, i.e., regime a in Fig. 3, shows
a similar scaling behavior as the flat geometry. At low activity this regime
spans over most of the time range until a steady state with low defect density
is reached, analogous to the behavior of the planar geometry (see the blue and
red lines in Fig. 5 [A] at $v_{0}=0.305$, where the black dot indicates the
transition towards the steady state configuration). The topological constraint
and the curvature of the sphere therefore do not affect the grain coarsening.
However, at high enough activity the dynamical process of coarsening changes
qualitatively, for which a transition to a faster coarsening regime becomes
evident. The transition into this regime happens earlier with increased
activity, as seen in Fig. 5 [A] ($v_{0}=0.5$, $0.425$, and $0.35$; blue
lines). Much faster coarsening processes, with much larger scaling exponents,
are observed, with a scaling behavior of roughly $t^{-2}$ (see Fig. 5 [A,C]).
These scaling exponents are calculated for the time range combining regimes b
and c of 3. The transition between these regimes strongly depends on the
specific configurations, as seen in 5 [B] showing results from all the
individual simulation runs on the sphere and the flat domains for $v_{0}=0.5$.
As the activity increases, the transition into the fast coarsening regime
occurs at earlier time, and the final configuration is reached earlier as well
(see the black diamonds in Fig. 5 [A] at $v_{0}=0.5$, $0.425$, and $0.35$,
which mark the start and end times for the fast coarsening regime, and the
black dots which mark the transition to the steady state).

###  3.5 Testing defect avalanche phenomenon

One could ask if the grain boundary collapse within this fast coarsening
regime shows any features of avalanche phenomenon. To explore this we follow
the characterization of avalanche event given in [77], through the number of
defects $N(t)$ and the defect velocity $v_{d}$, and a sudden change of both
quantities. Within an avalanche the defect velocity is significantly larger
than during the whole coarsening process and the number of defects reduces
faster. The later has already been identified by the time evolution of
$\mathcal{F}_{\mathbf{Q}}$ in Fig. 5, which relates to the defect density and
thus also to the number of defects. Here we examine the velocity of the
defects in this fast coarsening regime, with the comparison between the sphere
and plane domains for $v_{0}=0.5$.

Figure 6: Statistics of defect velocities: for sphere (top) and plane (bottom)
domains at a large activity $v_{0}=0.5$. [A] Evolution of average defect
velocity magnitude determined at each simulation timestep
$|v|(t)=\frac{1}{N(t)}\,\sum_{n=1}^{N(t)}\|v_{d}(t,n)\|$. This data
corresponds to Simulation 4 in Fig. 3 and is shown together with the
distortion energy $\mathcal{F}_{\mathbf{Q}}$ (blue or red curve, where
diamonds mark the times of snapshots shown in Fig. 2). [B] Ensemble average of
defect velocity magnitude
$[v](t)=\frac{1}{M}\sum_{m=1}^{M}\frac{1}{N_{m}(t)}\sum_{n=1}^{N_{m}(t)}\|v_{d}(t,n,m)\|$
across all simulations ($M=50$ for sphere geometry, $M=10$ for flat geometry),
together with the individual plots of $|v|(t)$ for all $M$ runs (blue and red
lines). [C] Probability distribution for defect velocity magnitude for
combined data set across all the simulations and the complete time domain
$t\in[10^{2},10^{5}]$ (black line), and for temporal binned data (blues and
red lines) where $[10^{2},10^{5}]$ is separated into $50$ slices of
logarithmic scaled boxes (as in [78]). The power law fitting, with the
obtained exponent $s$, for high velocities $|v_{d}|/\langle|v_{d}|\rangle\geq
1$ is marked by black dashed line.

Fig. 6 shows the evolution of the average defect velocity magnitude for
Simulation 4 of Fig. 3, its ensemble average over all simulations, and the
probability distribution of this data over the complete time domain and of
temporal binned data. All results are shown for the sphere and the plane
domain. For high velocities $|v_{d}|/\langle|v_{d}|\rangle\geq 1$ a power law
scaling with an exponent of roughly $-5.5$ can be identified for the
probability distribution for both the sphere and plane domains; see Fig. 6
[C]. This strong decay cannot be exclusively attributed to the fast coarsening
regime. First, the distributions are similar for the sphere and plane domains
at high velocities while the fast coarsening regime is not present for the
planar domain. Second, the probability distribution of the considered time
bins, i.e., the blue lines in Fig. 6 [C], shows essentially the same behaviour
as the overall distribution, which indicates that the high velocities are not
limited to any specific time span within the coarsening process. This is
confirmed by the time evolution of the averaged defect velocity magnitude in
Fig. 6 [A,B], which does not show a significantly larger velocity around the
time instance $t_{4}$ (fast coarsening regime). We only obtain a plateau of
lower velocities before the fast coarsening regime. Thus, these results are
not sufficient to identify the grain boundary collapse observed in the fast
coarsening regime as an avalanche phenomenon.

##  4 Discussion and Conclusions

The complex interplay of orientational and positional defects under
topological constraints, geometric effects, and active driving has been
investigated for active smectics on a sphere using an active surface PFC
model. Above a critical activity threshold we have identified a stable state
of rotating spirals which emerge from two rotating defect clusters. While the
number of defects in these clusters varies, they are always localized and
positioned to maximize the distance between each other. Such spiral states are
reminiscent to spirals in reaction-diffusions on a sphere [79, 80, 81, 82, 83,
84, 85].

The coarsening process that the system undergoes before it reaches this stable
state significantly differs from self-similar coarsening in passive smectics
or active smectics in the planar geometry. While several features of grain
boundary motion at early coarsening are similar to the passive case,
differences associated with activity are already present in this early regime.
These differences are rotating defect clusters, which lead to an enhanced
coarsening with defect density scaling of roughly $t^{-1/2}$ instead of
roughly $t^{-1/3}$ in the passive case. This increase in the coarsening rate
can be associated with activity and is independent of the considered geometry
of sphere or planar domain. A key characteristic for active smectics on a
sphere is the appearance of a fast coarsening regime. It is present at late
coarsening, just before the steady state is reached. This regime is
characterized by a scaling of roughly $t^{-2}$ for the time evolution of
defect density and is only present on the sphere domain. It is associated with
the collapse of grain boundary of a rotating grain. This collapse forms a new
defect cluster which is located on a pole position opposite to the already
present rotating defect cluster and rotates in the opposite direction, forming
the rotating spiral state. This process of grain boundary collapse and fast
coarsening is found to be different from the avalanche phenomenon. The
velocity of the defects within this regime is not significantly larger than
that during the whole coarsening process. Our simulations also show that this
fast coarsening regime follows a time period of frustration, associated with a
plateau or even an increase in the transient of the distortion energy
$\mathcal{F}_{\mathbf{Q}}$ and significantly lower defect velocities. The
extension of this plateau strongly depends on the local configuration.

This investigation on a sphere can only be the first step to understand the
complex interplay in the coarsening process of active surface smectics. For
surfaces with the same topology but varying curvature we speculate that a
similar steady state configuration will be reached. However, varying surface
curvature potentially localizes the rotating defect clusters at extrema of
surface curvature, e.g., in prolate configurations, or even leads to more
complicated metastable states for more complex geometries. We further
speculate that a similar collapse of defects will be also observed on more
general geometries just before the steady state configuration is reached.
However, the observed plateau in the transient of the distortion energy before
this event will probably strongly depend on local geometric features. All
these speculations ask for further numerical investigation and experimental
realization. Numerical methods which deal with non-spherical geometries and
surface PFC or surface liquid crystals (polar and nematic) have been developed
[86, 87, 88, 41] and need to be extended to the active case. Also, the
evolution of the surface in response to the surface liquid crystals on it is
an extension that can be studied numerically [42]. Extending this to surface
smectics will allow us to consider surface instabilities, such as the
Helfrich-Hurault elastic instability, which if considered for smectic shells
can lead to surface undulations and buckling [89]. Combined with activity this
might be the path to model and simulate some biological applications (e.g.,
thosed mentioned in [22]) as active surface smectics.

##  5 Methods

###  5.1 Numerical solution by (vector) spherical harmonics

We expand $\psi$ and $\mathbf{p}$ in (vector) spherical harmonics so that Eqs.
(5) and (6) reduce to a set of ordinary differential equations for the time-
dependent expansion coefficients of $\psi$ and $\mathbf{p}$. Let
$\mathcal{I}_{n}=\\{(l,m)\,:\,0\leq l\leq n,\,|m|\leq l\\}$ be an index set of
the spherical harmonics $Y^{m}_{l}:{\mathcal{S}}\to\mathbb{C}$ up to order
$n$. The scalar field $\psi$ can be expanded as

$\psi(\mathbf{r},t)=\\!\\!\\!\sum_{(l,m)\in\,\mathcal{I}_{\infty}}\\!\\!\\!\hat{\psi}_{lm}(t)Y_{l}^{m}(\mathbf{r}),$
(11)

with the expansion coefficients $\hat{\psi}_{lm}(t)$. For the vector field
$\mathbf{p}:{\mathcal{S}}\to\mathrm{T}{\mathcal{S}}$ we consider the
decomposition
$\mathbf{p}(\mathbf{r},t)={\operatorname{grad}_{{\mathcal{S}}}}p_{1}(\mathbf{r},t)+{\operatorname{rot}_{{\mathcal{S}}}}p_{2}(\mathbf{r},t)$,
and construct a tangent vector field basis from the gradient
${\operatorname{grad}_{{\mathcal{S}}}}$ and curl
${\operatorname{rot}_{{\mathcal{S}}}}=\mathbf{\hat{u}}\times{\operatorname{grad}_{{\mathcal{S}}}}$
of the spherical harmonics basis functions. We consider the vector spherical
harmonics
$\mathbf{y}_{lm}^{(1)}(\mathbf{r})=R\,{\operatorname{grad}_{{\mathcal{S}}}}Y_{l}^{m}(\mathbf{r})$
and
$\mathbf{y}_{lm}^{(2)}(\mathbf{r})=-\frac{\mathbf{r}}{\lVert\mathbf{r}\rVert}\times\mathbf{y}_{lm}^{(1)}(\mathbf{r})$,
and represent $\mathbf{p}$ by

$\mathbf{p}(\mathbf{r},t)=\sum_{i=1}^{2}\sum_{(l,m)\in\,\mathcal{I}_{\infty}}\\!\\!\\!\hat{p}_{lm}^{(i)}(t)\mathbf{y}_{lm}^{(i)}(\mathbf{r}),$
(12)

with expansion coefficients $\hat{p}_{lm}^{(i)}(t)$. The resulting Galerkin
scheme [90] reads

$\displaystyle\begin{split}\partial_{t}\hat{\psi}_{lm}(t)+\frac{l(l+1)}{R^{2}}\Big{\\{}{\varepsilon}+\Big{[}1-\frac{l(l+1)}{R^{2}}\Big{]}^{2}\Big{\\}}\hat{\psi}_{lm}(t)+\frac{l(l+1)}{R^{2}}\hat{\nu}_{lm}(t)-v_{0}\frac{l(l+1)}{R}\hat{p}_{lm}^{(1)}(t)&=0\,,\end{split}$
(13)
$\displaystyle\begin{split}\partial_{t}\hat{p}_{lm}^{(i)}(t)+\Big{[}\frac{l(l+1)}{R^{2}}+D_{r}\Big{]}\big{(}C_{1}\hat{p}_{lm}^{(i)}(t)+C_{2}\hat{q}_{lm}^{(i)}(t)\big{)}+v_{0}\frac{\delta_{i1}}{R}\hat{\psi}_{lm}(t)&=0,\end{split}$
(14)

with $i\in\\{1,2\\}$, $(l,m)\in\mathcal{I}_{n}$, and
$t\in[t_{0},t_{\mathrm{end}}]$. Here $\hat{\nu}_{lm}(t)$ are the expansion
coefficients of $\nu=\psi^{3}$, $\hat{q}_{lm}^{(i)}(t)$ are the expansion
coefficients of $\mathbf{q}=\lVert\mathbf{p}\rVert^{2}\mathbf{p}$, and
$t_{\mathrm{end}}$ represents the simulated time range starting at $t_{0}=0$.
Quadrature on the sphere ${\mathcal{S}}$ is realized by evaluating $\psi$ and
$\mathbf{p}$ in Gaussian points $\\{(\theta_{i},\phi_{j})\,:\,1\leq i\leq
N_{\theta},1\leq j\leq N_{\phi}\\}$, where $N_{\theta}$ and $N_{\phi}$ are the
numbers of grid points along the polar and azimuthal coordinates,
respectively. An appropriate quadrature rule [91] is applied. We use a second-
order accurate scheme similar to that described in [92] for the time-
discretization of Eqs. (13) and (14). The implementation of the vector
spherical harmonics is based on the toolbox SHTns [91] and identical to the
one used in [23].

For numerical solution procedure we discretize the spherical domain with
$N_{\theta}=512$ and $N_{\phi}=1024$, and used spherical harmonics expansion
of order $n=500$ and a time step $\tau=0.05$ for time-discretization. The
simulations are performed on the time domain $[0,10^{5}]$.

###  5.2 Identification and tracking of disclinations and dislocations

To detect disclinations and dislocations in the smectic pattern we use the
nematic order parameter of Eq. (9) and its distortion energy defined in Eq.
(10). For these quantities several techniques have been proposed (see, e.g.,
[93]) to identify dislocations and disclinations. Here we use the approach to
identify defects as localized maxima of distortion energy in Eq. (10). It is
noted that in this context the distinction between dislocation and
disclination is dropped. To exclude faulty detections due to numeric
fluctuations in $f_{d}(\mathbf{Q})$ we use a lower threshold
$f_{d}(\mathbf{Q})>1$. Furthermore we filter the resulting positions by a
distance threshold of
$\|\mathbf{x}_{i}-\mathbf{x}_{j}\|_{\mathbb{R}^{3}}>\pi/2$ to reduce possible
sets of maxima in $f_{d}(\mathbf{Q})$ close to a disclination to a single
position.

At the next step we aim to identify larger coherent structures, like rotating
defect clusters and grain boundaries. We evaluate a hierarchical clustering
[94] of the positions along Euclidean distance with the requirement of single
connection, using the implementation provided by scipy.cluster [95]. The
coherent structures are given by a slice (of the defect position dendrogram)
with threshold $3/2\pi$. As the final step we connect the instantaneous defect
positions and clusters by using the trackpy [96] particle tracking library.
From these trajectories we obtain defect velocity magnitude by central
difference
$v_{d}(t^{i})=[\|\mathbf{x}^{i}-\mathbf{x}^{i-1}\|/(t^{i}-t^{i-1})+\|\mathbf{x}^{i+1}-\mathbf{x}^{i}\|/(t^{i+1}-t^{i})]/2$.

###  5.3 Fast coarsening regime and steady state regime tracking

To quantify and track the fast coarsening regime we consider the distortion
energy transient $(t,\mathcal{F}_{\mathbf{Q}}(t))$ in the $\log-\log$ space
$(\underline{t}=\log_{10}t$,
$\underline{\mathcal{F}}_{\mathbf{Q}}=\log_{10}\mathcal{F}_{\mathbf{Q}})$.
There we separate the $\log-\log$ transient in equidistant slices
$[\underline{t}_{i},\underline{t}_{i}+\delta t],\,i\in\\{0,\ldots 250\\}$,
evaluate the relative energy decrease in each slice
$d_{i}=[\underline{\mathcal{F}}_{\mathbf{Q}}(\underline{t}_{i}+\delta
t)-\underline{\mathcal{F}}_{\mathbf{Q}}(\underline{t}_{i})]/\max_{k}[\underline{\mathcal{F}}_{\mathbf{Q}}(\underline{t}_{k}+\delta
t)-\underline{\mathcal{F}}_{\mathbf{Q}}(\underline{t}_{k})]$ and a threshold
by $d_{i}>0.75$ to obtain time slices containing fast coarsening. In a second
step neighbouring slices with fast coarsening regime are combined.

Applying this method to the transients of averaged energies, as described in
Fig. 5, we obtain a single time domain with fast coarsening regime. Note that
when applying such a method to an energy transient of a single simulation run
we have to reduce noise in the transient, e.g., by convolution with a Gaussian
function along temporal axis in the $\log-\log$ space, and then yield possibly
more than one time domains with fast coarsening regime.

To identify a steady state regime we determine the average nematic distortion
energy $\langle\mathcal{F}_{\mathbf{Q}}\rangle$ in the last $5\%$ of
simulation time domain $[0,T]$. If the energy at $t=0.95T$, the begin of
averaging domain, is less than a value $5\%$ higher than the averaged energy,
i.e., if
$\mathcal{F}_{\mathbf{Q}}(0.95T)<1.05\langle\mathcal{F}_{\mathbf{Q}}\rangle$,
we classify the simulation as reaching the steady state. If a steady state is
reached, we find the steady state domain $[T_{s},T]$ by traversing the
simulation time steps backward, starting from $t=0.95T$, as long as the
criteria
$\mathcal{F}_{\mathbf{Q}}(t)<1.05\langle\mathcal{F}_{\mathbf{Q}}\rangle$ is
satisfied.

###  5.4 Decay rate fitting

To estimate the power-law decay rates as done in Figures 5 and 6, we transform
the considered data into the $\log-\log$ space and perform linear regression
on the relevant subset of the data to obtain the decay rate $r$. For Fig. 5 we
use the full data set in cases where no steady state was tracked, while in
cases with steady state we consider data sets either restricted to $[0,T_{s}]$
(if no fast coarsening regime was tracked) or restricted to the regime of fast
coarsening. In Fig. 6 we restrict the data set to velocities faster than the
average, i.e., $|v_{d}|/\langle|v_{d}|\rangle\geq 1$.

###  5.5 Statistics of defect velocities

As we are concerned with evolution in a time span across several orders of
magnitude, we typically describe the phenomena in $\log-\log$ space. To reduce
noise in the simulation results we apply the temporal binning with $\log$
sized bins as provided by numpy.logspace[97] such that the considered bins are
equally sized in $\log$ space.

To characterize the defect velocity we apply several notions reflecting the
possible perspectives on the data. As the first notion we consider

$|v|(t)=\frac{1}{N(t)}\sum_{n=1}^{N(t)}\|v_{d}(t;n,I)\|,$

the average defect velocity magnitude across a defect configuration
(consisting of $N(t)$ defects) at a certain time $t$ and in a simulation $I$.
To obtain a characteristic defect velocity magnitude independent of
instantaneous defect configuration and specific realization of evolution,
which highly depend on the random initial value, we introduce an ensemble
average across $M$ simulations to obtain the characteristic average defect
velocity magnitude

$[v](t)=\frac{1}{M}\sum_{m=1}^{M}\frac{1}{N_{m}(t)}\sum_{n=1}^{N_{m}(t)}\|v_{d}(t,n,m)\|.$

Furthermore, we consider the probability distribution of defect velocities to
obtain a notion compatible to [78]. For this purpose we consider the magnitude
of defect velocity as a random variable. We calculate the associated expected
value $\langle|v_{d}|\rangle$ and probability distribution numerically from
the combined data set of all defect configurations across all simulations,
$\|v_{d}(t,n,m)\|,\forall n=1\ldots N,\,m=1\ldots M,\,t\in[10^{2},10^{5}]$, as
shown in Fig. 6 [C] (black curves). To test whether the observed different
coarsening regimes have an impact on this probability distribution we repeat
the calculation by restricting the defect velocity data to different temporal
bins with the use of 50 slices of logarithmic scaled boxes.

## References

## References

  * [1] de Gennes P G and Prost J 1993 The Physics of Liquid Crystals (Second Edition, Oxford Science Publications)
  * [2] Chaikin P M and Lubensky T C 1995 Principles of condensed matter physics (Cambridge University Press)
  * [3] Elder K R, Vinals J and Grant M 1992 Phys. Rev. Lett. 68 3024–3027
  * [4] Elder K R, Vinals J and Grant M 1992 Phys. Rev. A 46 7618–7629
  * [5] Cross M C and Meiron D I 1995 Phys. Rev. Lett. 75 2152–2155
  * [6] Boyer D and Vinals J 2001 Phys. Rev. E 64 050101
  * [7] Boyer D and Vinals J 2002 Phys. Rev. E 65 046119
  * [8] Abukhdeir N M and Rey A D 2008 New J. Phys. 10 063025
  * [9] Harrison C, Adamson D, Cheng Z, Sebastian J, Sethuraman S, Huse D, Register R and Chaikin P 2000 Science 290 1558–1560
  * [10] Purvis L and Dennin M 2001 Phys. Rev. Lett. 86 5898–5901
  * [11] Pershan P S and Prost J 1975 J. Appl. Phys. 46 2343–2353
  * [12] Harrison C, Cheng Z, Sethuraman S, Huse D, Chaikin P, Vega D, Sebastian J, Register R and Adamson D 2002 Phys. Rev. E 66 011706
  * [13] Lopez-Leon T, Fernandez-Nieves A, Nobili M and Blanc C 2011 Phys. Rev. Lett. 106 247802
  * [14] Liang H L, Noh J, Zentel R, Rudquist P and Lagerwall J P 2013 Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 371 20120258
  * [15] Bechinger C, Di Leonardo R, Löwen H, Reichhardt C, Volpe G and Volpe G 2016 Rev. Mod. Phys. 88 045006
  * [16] Adhyapak T C, Ramaswamy S and Toner J 2013 Phys. Rev. Lett. 110 118102
  * [17] Chen L and Toner J 2013 Phys. Rev. Lett. 111 088701
  * [18] Romanczuk P, Chate H, Chen L, Ngo S and Toner J 2016 New J. Phys. 18 063015
  * [19] Saha S, Agudo-Canalejo J and Golestanian R 2020 Phys. Rev. X 10 041009
  * [20] Tarama S, Egelhaaf S U and Löwen H 2019 Phys. Rev. E 100 022609
  * [21] Doostmohammadi A, Ignes-Mullol J, Yeomans J M and Sagues F 2018 Nat. Commun. 9 3246
  * [22] Jülicher F, Prost J and Toner J 2022 Phys. Rev. E 106 054607
  * [23] Praetorius S, Voigt A, Wittkowski R and Löwen H 2018 Phys. Rev. E 97 052615
  * [24] Huang Z F, Löwen H and A V 2022 Commun. Phys. 5 294
  * [25] Bowick M J and Giomi L 2009 Adv. Physics 58 449
  * [26] Thomson J 1904 Philos. Mag. 7 237–265
  * [27] Erber T and Hockney G 1991 J. Phys. A 24 L1369–L1377
  * [28] Smale S 1998 Math. Intelligencer 20 7–15
  * [29] Bowick M J, Nelson D R and Travesset A 2000 Phys. Rev. B 62 8738–8751
  * [30] Bausch A, Bowick M, Cacciuto A, Dinsmore A, Hsu M, Nelson D, Nikolaides M, Travesset A and Weitz D 2003 Science 299 1716–1718
  * [31] Backofen R, Voigt A and Witkowski T 2010 Phys. Rev. E 81 025701(R)
  * [32] Lubensky T and Prost J 1992 J. Phys. II 2 371–382
  * [33] Turner A M, Vitelli V and Nelson D R 2010 Rev. Mod. Phys. 82 1301–1348
  * [34] Dzubiella J, Schmidt M and Löwen H 2000 Phys. Rev. E 62 5081–5091
  * [35] Bates M A, Skacej G and Zannoni C 2010 Soft Matter 6 665
  * [36] Shin H, Bowick M J and Xing X 2008 Phys. Rev. Lett. 101 037802
  * [37] Dhakal S, Solis F J and Olvera de la Cruz M 2012 Phys. Rev. E. 86 011709
  * [38] Koning V, Lopez-Leon T, Fernandez-Nieves A and Vitelli V 2013 Soft Matter 9 4993
  * [39] Napoli G and Vergori L 2012 Phys. Rev. Lett. 108 207803
  * [40] Henkes S, Marchetti M C and Sknepnek R 2018 Phys. Rev. E 97 042605
  * [41] Nitschke I, Nestler M, Praetorius S, Löwen H and Voigt A 2018 Proc. Roy. Soc. A 474 20170686
  * [42] Nitschke I, Reuther S and Voigt A 2020 Proc. Roy. Soc. A 476 20200313
  * [43] Smallenburg F and Löwen H 2016 J. Chem. Phys. 144
  * [44] Mosna R A, Beller D A and Kamien R D 2012 Phys. Rev. E 86 011707
  * [45] Xing X, Shin H, Bowick M J, Yao Z, Jia L and Li M H 2012 Proc. Nat. Acad. Sci. (USA) 109 5202
  * [46] Allahyarov E, Voigt A and Löwen H 2017 Soft Matter 13 8120–8135
  * [47] Marchetti M C, Joanny J F, Ramaswamy S, Liverpool T B, Prost J, Rao M and Simha R A 2013 Rev. Mod. Phys. 85 1143
  * [48] Berry J, Brangwynne C P and Haataja M 2018 Rep. Prog. Phys. 81 046601
  * [49] Bowick M J, Fakhri N, Marchetti M C and Ramaswamy S 2022 Phys. Rev. X 12 010501
  * [50] Shankar S, Souslov A, Bowick M J, Marchetti M C and Vitelli V 2022 Nature Rev. Phy. 4 380–398
  * [51] Cross M C and Hohenberg P C 1993 Rev. Mod. Phys. 65 851–1112
  * [52] Keber F C, Loiseau E, Sanchez T, DeCamp S J, Giomi L, Bowick M J, Marchetti M C, Dogic Z and Bausch A R 2014 Science 345 1135–1139
  * [53] Alaimo F, Koehler C and Voigt A 2017 Sci. Rep. 7 5211
  * [54] Apaza L and Sandoval M 2018 Soft Matter 14 9928–9936
  * [55] Pearce D J G, Ellis P W, Fernandez-Nieves A and Giomi L 2019 Phys. Rev. Lett. 122 168002
  * [56] Nestler M and Voigt A 2022 Commun. Comput. Phys. 31 947–965
  * [57] Sknepnek R and Henkes S 2015 PHys. Rev. E 91 022306
  * [58] Nitschke I, Reuther S and Voigt A 2019 Phys. Rev. Fluids 4 044002
  * [59] Tanner K, Mori H, Mroue R, Bruni-Cardoso A and Bissell M J 2012 Proc. Nat. Acad. Sci. (USA) 109 1973–1978
  * [60] Wang H, Lacoche S, Huang L, Xue B and Muthuswamy S K 2013 Proc. Nat. Acad. Sci. (USA) 110 163–168
  * [61] Happel L, Wenzel D and Voigt A 2022 EPL 138 67002
  * [62] Tan T H, Amiri A, Seijo-Barandiarán I, Staddon M F, Materne A, Tomas S, Duclut C, Popović M, Grapin-Botton A and Jülicher F 2022 bioRxiv
  * [63] Achim C V, Wittkowski R and Löwen H 2011 Phys. Rev. E 83 061712
  * [64] Radtke L and Voorhees P W 2002 Growth and Coarsening: Ostwald Ripening in Material Processing (Springer)
  * [65] Backofen R, Gräf M, Potts D, Praetorius S, Voigt A and Witkowski T 2011 Multiscale Modeling & Simulation 9 314–334
  * [66] Backofen R, Barmak K, Elder K and Voigt A 2014 Acta Mater. 72–77
  * [67] Menzel A M and Löwen H 2013 Phys. Rev. Lett. 110 055702
  * [68] Menzel A M, Ohta T and Löwen H 2014 Phys. Rev. E 89 022301
  * [69] Alaimo F, Praetorius S and Voigt A 2016 New J. Phys. 18 083008
  * [70] Alaimo F and Voigt A 2018 Phys. Rev. E 98 032605
  * [71] Ophaus L, Gurevich S and Thiele U 2018 Phys. Rev. E 98 022608
  * [72] Huang Z F, Menzel A M and Löwen H 2020 Phys. Rev. Lett. 125 218002
  * [73] Ophaus L, Knobloch E, Gurevich S V and Thiele U 2021 Phys. Rev. E 103 032601
  * [74] Löwen H 2010 J. Phys.: Condensed Matter 22 364105
  * [75] Koehler C, Backofen R and Voigt A 2016 Phys. Rev. Lett. 116 135502
  * [76] Menzel A M 2015 Phys. Rep. 554 1–45
  * [77] Chan P Y, Tsekenis G, Dantzig J, Dahmen K A and Goldenfeld N 2010 Phys. Rev. Lett. 105 015502
  * [78] Angheluta L, Jeraldo P and Goldenfeld N 2012 Phys. Rev. E 85 011153
  * [79] Maselko J and Showalter K 1989 Nature 339 609–611
  * [80] McQuillan P and Gomatam J 1996 J. Phys. Chem. 100 5157–5159
  * [81] Gomatam J and Amdjadi F 1997 Phys. Rev. E 56 3913–3919
  * [82] Davidsen J, Glass L and Kapral R 2004 Phys. Rev. E 70 056203
  * [83] Amdjadi F and Gomatam J 2005 J. Comput. Appl. Math. 182 472–486
  * [84] Calhoun D A, Helzel C and LeVeque R J 2008 SIAM Review 50 723–752
  * [85] Sigrist R and Matthews P 2011 SIAM J. Appl. Dyn. Sys. 10 1177–1211
  * [86] Aland S, Lowengrub J and Voigt A 2011 Phys. Fluids 23 062103
  * [87] Aland S, Lowengrub J and Voigt A 2012 Phys. Rev. E 86 046321
  * [88] Nestler M, Nitschke I, Praetorius S and Voigt A 2018 Journal of Nonlinear Science 28 147–191
  * [89] Blanc C, Durey G, Kamien R D, Lopez-Leon T, Lavrentovich M O and Tran L 2023 Rev. Mod. Phys. 95 015004
  * [90] Hesthaven J S and Warburton T 2007 Nodal discontinuous Galerkin methods: algorithms, analysis, and applications (Springer Science & Business Media)
  * [91] Schaeffer N 2013 Geochemistry, Geophysics, Geosystems 14 751–758
  * [92] Backofen R, Gräf M, Potts D, Praetorius S, Voigt A and Witkowski T 2011 Multiscale Modeling & Simulation 9 314–334
  * [93] Wenzel D, Nestler M, Reuther S, Simon M and Voigt A 2021 Computational Methods in Applied Mathematics
  * [94] Müllner D et al. 2013 Journal of Statistical Software 53 1–18
  * [95] Virtanen P, Gommers R, Oliphant T E, Haberland M, Reddy T, Cournapeau D, Burovski E, Peterson P, Weckesser W, Bright J, van der Walt S J, Brett M, Wilson J, Millman K J, Mayorov N, Nelson A R J, Jones E, Kern R, Larson E, Carey C J, Polat İ, Feng Y, Moore E W, VanderPlas J, Laxalde D, Perktold J, Cimrman R, Henriksen I, Quintero E A, Harris C R, Archibald A M, Ribeiro A H, Pedregosa F, van Mulbregt P and SciPy 10 Contributors 2020 Nature Methods 17 261–272
  * [96] Allan D, van der Wel C, Keim N, Caswell T, Wieker D, Verweij R, Reid C, Grueter L, Ramos K, Perry R et al. 2019 soft-matter/trackpy: Trackpy v0. 4.2
  * [97] Harris C R, Millman K J, van der Walt S J, Gommers R, Virtanen P, Cournapeau D, Wieser E, Taylor J, Berg S, Smith N J, Kern R, Picus M, Hoyer S, van Kerkwijk M H, Brett M, Haldane A, del Río J F, Wiebe M, Peterson P, Gérard-Marchant P, Sheppard K, Reddy T, Weckesser W, Abbasi H, Gohlke C and Oliphant T E 2020 Nature 585 357–362

## Acknowledgements

We are thankful for the support of this work by the Deutsche
Forschungsgemeinschaft (DFG) through the grants LO 418/20-2 (for H.L.) and VO
899/19-2 (for A.V.), PR 1705/1-1 (for S.P.) within the DFG Research Unit 3013,
and the National Science Foundation (NSF) under Grant No. DMR-2006446 (for
Z.-F.H.). We further acknowledge computing resources provided at FZ Jülich
within grant pfamdis and at ZIH/TU Dresden within grant WIR.

## Author contributions statement

Z.-F.H., H.L. and A.V. developed the model and conceived the computational
experiments. S.P. implemented the simulation algorithm for the spherical
geometry, and Z.-F.H implemented the simulation algorithm for the planar case.
M.N. implemented the postprocessing algorithms and conducted the computational
experiments. M.N. and A.V. analysed the results and wrote the manuscript. All
authors reviewed and revised the manuscript.

## Additional information

Accession codes: Data is available upon reasonable request.
Competing interests: There are no competing interests.

## SI

### Animations of Evolution

To supplement the description for the evolution of smectic defect density
distributions presented in Figures 2 and 4, we provide visualization
animations which can be found at

https://datashare.tu-dresden.de/s/jXJ5yBNg3wKZF94

#### EvolSmecticLines

This video corresponds to 2, showing the rescaled density field $\psi$ in the
front and back perspectives on a linear time scale.

#### D450_V0.5

Here we combine the defect distributions along 100 logarithmic temporal bins
to an animation. On the left side the transient of nematic distortion energy
is depicted, while the red diamonds indicate the domain of the temporal bin.
The right side provides the defect density distribution for the time bin.
These results correspond to the snapshots presented in 4 [A].

#### R100_V0.5

Analog to D450_V0.5, this animation corresponds to the result presented in 4
[B]. The defect densities have been derived from simulation results shown in
EvolSmecticLines.