Active smectics on a sphere

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.


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 α for the characteristic domain size with time t and a scaling exponent α.The exponent α 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 [8,12] 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 [7][8][9]12]).
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.

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.

Active liquid crystals on a sphere
The active liquid crystalline phases, which have been of tremendous interest in recent years [15,[47][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.

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 [24,[67][68][69][70][71][72][73] have been introduced.These models are related to the PFC model for liquid crystals [63] and have been formulated on a sphere in [23,74] to study active crystals.

Results
In the following we first introduce the surface active PFC model and then describe largescale 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.

Model
In [23] the active PFC model introduced in [67] has been formulated on a sphere S = R S 2 with radius R, where S 2 is the two dimensional unit sphere embedded in 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 r on the sphere S is parametrized by r(θ, φ) = Rû(θ, φ) with the orientational unit vector û(θ, φ) = (sin θ cos φ, sin θ sin φ, cos θ) T and the spherical coordinates θ ∈ [0, π] and φ ∈ [0, 2π).We consider a rescaled density field ψ(r, t) and a polarization field p(r, t), which is tangential to S at r, i.e., p(r, t) = p θ (r, t)∂ θ û + p φ (r, t)∂ φ û ∈ T r S with scalar polarization fields p θ (r, t) and p φ (r, t) and the tangent space T r S of the sphere S at point r.
The free-energy functionals read with ǫ < 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., and the surface Laplace-Beltrami operator is defined as △ S = div S grad S .The dynamic equations describe the active-particle transport tangential to S and read with v 0 the strength of particle self-propulsion and D r the rotational diffusion constant.The vector-Laplacian △ dR is the surface Laplace-de Rham operator △ dR = −grad S div S − rot S Rot S , where 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 → ∞.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 ǫ = −0.98 and the average density ψ 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 [23,76] 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 ∈ {0.305, 0.35, 0.425, 0.5}.Each simulation starts from uniform ψ(0) = ψ 0 disturbed by random noise ψ ′ with magnitude |ψ ′ | < 0.01 applied at each grid point.
In order to identify orientational defects we define a nematic order parameter almost everywhere x : grad S ψ(x) = 0: with Π(x) = I − ν(x)ν(x) the surface identity of the considered geometry and ν 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 Q w.r.t the coordinates of embedding space R 3 which will be used to identify defects.For details we refer to the methods section 5.

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 ψ.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. 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 F Q , together with the considered time instances, is shown in Fig. 3 (top panel of Simulation 4).

Distinct coarsening regimes
The transients of the distortion energy F 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 , . . ., 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×100 spatial bins) of the positions.The resulting twodimensional 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 θ, x 1 to azimuthal angle φ) 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 F 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.

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 F Q over time, which can be related to the statistics on the time decay of defect density.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).

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 F 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.
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 |/ |v d | ≥ 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.

Discussion and Conclusions
The complex interplay of orientational and positional defects under topological constraints, geometric effects, and active driving has been investigated for active This data corresponds to Simulation 4 in Fig. 3 and is shown together with the distortion energy F Q (blue or red curve, where diamonds mark the times of snapshots shown in Fig. 2).[B] Ensemble average of defect velocity magnitude 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 F 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 nonspherical geometries and surface PFC or surface liquid crystals (polar and nematic) have been developed [41,[86][87][88] 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.

Numerical solution by (vector) spherical harmonics
We expand ψ and 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 ψ and p.Let I n = {(l, m) : 0 ≤ l ≤ n, |m| ≤ l} be an index set of the spherical harmonics Y m l : S → C up to order n.The scalar field ψ can be expanded as with the expansion coefficients ψlm (t).For the vector field p : S → TS we consider the decomposition p(r, t) = grad S p 1 (r, t) + rot S p 2 (r, t), and construct a tangent vector field basis from the gradient grad S and curl rot S = û × grad S of the spherical harmonics basis functions.We consider the vector spherical harmonics y with expansion coefficients p(i) lm (t).The resulting Galerkin scheme [90] reads with i ∈ {1, 2}, (l, m) ∈ I n , and t ∈ [t 0 , t end ].Here νlm (t) are the expansion coefficients of ν = ψ 3 , q(i) lm (t) are the expansion coefficients of q = p 2 p, and t end represents the simulated time range starting at t 0 = 0. Quadrature on the sphere S is realized by evaluating ψ and p in Gaussian points {(θ i , φ j ) : 1 ≤ i ≤ N θ , 1 ≤ j ≤ N φ }, where N θ and N φ 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 θ = 512 and N φ = 1024, and used spherical harmonics expansion of order n = 500 and a time step τ = 0.05 for time-discretization.The simulations are performed on the time domain [0, 10 5 ].

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 (Q) we use a lower threshold f d (Q) > 1. Furthermore we filter the resulting positions by a distance threshold of x i − x j R 3 > π/2 to reduce possible sets of maxima in f d (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π.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

Fast coarsening regime and steady state regime tracking
To quantify and track the fast coarsening regime we consider the distortion energy transient (t, F Q (t)) in the log − log space (t = log 10 t, F Q = log 10 F Q ).There we separate the log − log transient in equidistant slices [t i , t i + δt], i ∈ {0, . . .250}, evaluate the relative energy decrease in each slice 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 F Q 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 F Q (0.95T ) < 1.05 F Q , 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 F Q (t) < 1.05 F Q is satisfied.

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 |/ |v d | ≥ 1.

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 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 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 |v d | and probability distribution numerically from the combined data set of all defect configurations across all simulations, v d (t, n, m) , ∀n = 1 . . .N, m = 1 . . .M, t ∈ [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.
EvolSmecticLines This video corresponds to 2, showing the rescaled density field ψ in the front and back perspectives on a linear time scale.

Fig. 1
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 .
Figure 1.Smectic pattern formation and dynamics in sphere domain.[A] Snapshots of ψ 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 ψ with time steps of △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.

Figure 2 .
Figure 2. Grain growth during the coarsening process for v 0 = 0.5 in sphere domain, showing as a sequence of snapshots depicting ψ at various times marked in Simulation 4 of Fig. 3, with front view (top panels) and back view (bottom) in terms of 180 • 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.

Figure 3 .
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 F 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 .
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 , • • • , t 5 correspond to the transients marked in Fig. 3 (Simulation 4) showing the nematic distortion energy F 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.

Figure 5 .
Figure 5. Average behavior of coarsening for sphere and plane domains.[A] The nematic distortion energy F 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 F 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.

Figure 6 .
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) = 1 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 ∈ [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 |/ |v d | ≥ 1 is marked by black dashed line.

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.