# Self-organization and shape change by active polarization in nematic droplets

Fabian Jan Schwarzendahl, Pierre Ronceray, Kimberly Weirich, Kinjal Dasbiswas

► To cite this version:

Fabian Jan Schwarzendahl, Pierre Ronceray, Kimberly Weirich, Kinjal Dasbiswas. Self-organization and shape change by active polarization in nematic droplets. Physical Review Research, 2021, 3 (4), pp.043061. 10.1103/PhysRevResearch.3.043061. hal-03623792

HAL Id: hal-03623792

https://hal.science/hal-03623792v1

Submitted on 2 Sep 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Self-organization and shape change by active polarization in nematic droplets

Fabian Jan Schwarzendahl,¹,² Pierre Ronceray,³ Kimberly L. Weirich,⁴ and Kinjal Dasbiswas¹

¹Department of Physics, University of California Merced, Merced, California 95343, USA

²Institut für Theoretische Physik II: Weiche Materie,
Heinrich-Heine-Universität Düsseldorf, 40225 Düsseldorf, Germany

³Aix Marseille Univ, Université de Toulon, CNRS, CPT,
Turing Center for Living Systems, Marseille, France

⁴Department of Materials Science and Engineering, Clemson University, Clemson, SC 29634, USA

(Dated: February 16, 2021)

Active forces occurring within cells can drive crucial biological processes that involve spontaneous organization and shape change, such as cell division. Motivated by recent *in vitro* experiments of nematic droplets of cytoskeletal filaments and motors that self-organize and divide, we present a minimal hydrodynamic model that combines the nonequilibrium kinetics of motor-filament interactions with equilibrium nematic phase separation. The motors organize within droplets and structure filaments into polarized aster defects. At large motor activity, they can even deform or divide the droplet, or form multi-aster chains of droplets. Our predicted phase diagram recapitulates these experimentally observed shapes.

***Introduction.*** Active mechanical forces enable living systems, particularly animal cells, to move, change shape, organize components, and divide. Subcellular cytoskeletal assemblies, comprising polar filaments and molecular motors that transduce biochemical reactions to generate active mechanical forces, drive these processes [1]. Understanding the general physical principles of living matter provides insight into cell biology, as well as guides the engineering of artificial cells that exhibit spatiotemporal organization of components and spontaneous shape change characteristic of cell division [2]. Model *in vitro* systems of purified cytoskeletal proteins, which capture elements of cell biological phenomena with only a fraction of the biochemical complexity occurring *in vivo*, exhibit a rich array of collective phenomena [3, 4] that motivates bio-inspired active matter theory [5].

Recently, phase segregated macromolecular droplets have emerged as model systems to investigate spatiotemporal organization in biological cells and phase separation has been proposed as a possible primitive means of subcellular organization in protocells [6]. These macromolecular liquids are typically composed of disordered proteins and nucleic acids, and consequently, internal droplet order as well as motor activity are absent. In contrast, recent experiments indicate that biopolymers with high aspect ratio, such as cytoskeletal filaments, can form phase separated droplets with orientational order, as in liquid crystals, because of the alignment of filaments in the dense phase [7]. This nematic order confers an equilibrium spindle shape to these droplets [8, 9], known as tactoids, which arises from a competition of droplet surface tension, the tendency of the filaments to align with the interface, and elasticity arising from bulk nematic order [10, 11].

The ordered structure in biomolecular fluids influences the emergence of collective phenomena in active systems [5]. When confined to a droplet, active forces lead to non-equilibrium phenomena such as droplet shape change [12–19], motility [20, 21] and dynamics governed

by the geometry of the confining droplet [22]. In addition to these active nematic fluids, the directed “walking” of motors on filaments, along with motor-based filament crosslinking, can lead to polar order, where filaments prefer to point in the same direction at the mesoscopic scale, and form self-organized defect structures, such as asters and vortices [23]. Such polar active states are fruitfully described by hydrodynamic theories [24–30]. In recent experiments, myosin motors were shown to self-organize at the midplane of the aforementioned actin-filament based nematic droplets [31]. When sufficiently active myosin motors are present, they deform the droplet, even splitting it into two. A simple free energy-based model of a nematic droplet, considering only the mutual alignment of the filaments and motors and adhesion of the droplet to the motor complex, was invoked to capture these key behaviors, but relied on arguments specific to the shape and structure of both the motor complexes and the droplet [8]. On the other hand, active mechanical forces, specifically the directed sliding of filaments by motors leading to their sorting by polarity [1], are expected to play a role in the dynamics of the organization. The effects of such active emergence of local polar order and associated defects, within a nematic droplet with a preferred orientational axis, have not yet been theoretically explored.

In this Letter, we combine continuum modeling that describes the structure of equilibrium nematic droplets with an active mechanical model for how motors move and slide filaments according to polarity, while self-organizing into localized asters. Using numerical simulations, complemented by theoretical analysis, we show that the resulting motor-filament self-organization destabilize droplets, giving rise to a rich array of experimentally observed structures including deformed, divided and multi-lobed droplets that can be generated by tuning one of two motor activity-dependent parameters.

***1D polarity sorting model.*** To build up intuition on the role of active forces in motor self-organization, we first
examine a simplified one dimensional setup of motors interacting with filaments, as sketched in Fig. 1(i). Here, red lines depict actin filaments, blue circles show myosin II motors, while the black arrows indicate the motors respective direction of motion. This scenario arises, for instance, in contractile actin bundles [32, 33] and when actin filaments are locally oriented along a nematic director [34]. Since motors walk towards a specific end of the polar filament (the barbed or plus end), we consider a density of left and right pointing filaments, denoted by $n_+$ and $n_-$ respectively. Momentum balance for the motor-filament system requires that a right (left) pointing filament is pushed to the right (left) by the motor as the motor moves in a single direction along the polar filament. The active motion of motors of density $m$ thus gives rise to mass fluxes of both motors and filaments. In the dilute regime, the fluxes of right and left pointing filaments can be expressed as $J_+ = \zeta n_+ m$ and $J_- = -\zeta n_- m$ respectively, where $\zeta$ is a parameter related to the active motor force. The net flux of motors is written as $-v_0(n_+-n_-)m$, where $v_0$ is the self-propulsion velocity of motors.

Including the diffusion of filaments and motors with coefficients $D$ and $D_m$ respectively, the flux conservation equations can be written in terms of filament density $\rho = n_+ + n_-$ and the 1D polarization $p = (n_+ - n_-)/\rho_0$ of filaments, where $\rho_0 = \langle \rho \rangle$ is the average density, giving

$$ \partial_t \rho = D \partial_x^2 \rho - \zeta \rho_0 \partial_x (pm), \quad (1) $$

$$ \partial_t p = D \partial_x^2 p - \zeta / \rho_0 \partial_x (\rho m), \quad (2) $$

$$ \partial_t m = D_m \partial_x^2 m + v_0 \partial_x (mp). \quad (3) $$

<figure><img src="image_2.png" /><figcaption>FIG. 1. Illustration of steps in the active polarity-sorting model. (i) Motors “walk” towards barbed ends of filaments which they slide in the opposite direction as a result of momentum conservation. (ii) Polarization and (iii) motor density at steady state as a function of position in a one-dimensional model (Eqs. S7-S6) that shows motor localization to center of bundle by polarity sorting. (iv) In the 2D nematic droplet model, the filaments are initially oriented along the long axis without any polarity preference. Motors bind to filaments and advect them according to their polarity. (v) At steady state, motors gather at the droplet midplane and sort filaments into an aster that pinches the droplet.</figcaption></figure>

The steady-state solution of Eqs. (2)-(3), (see SI I-II), is the one dimensional equivalent of an aster, where the populations of right and left pointing filaments are completely sorted to the right and the left of the aster, and is shown in Fig. 1(ii)-(iii). Note that in this minimal model, we have not considered the usual role of myosin motors as active crosslinkers which create pairwise forces on anti-parallel filaments, leading to self-straining flows proportional to local polarization [35]. This simple model reproduces the key experimental observation that motors migrate to the center of the droplet, and predicts that they induce strong polarity sorting.

**2D nematic droplet with motor-induced polarization.** To explore this prediction of active polarity sorting in nematic droplets and its implications for droplet shape, we build a more realistic hydrodynamic description for a suspension of filaments and motors on a frictional substrate that damps out large-scale fluid flows. In contrast with a thin fluid film, the filaments we seek to describe aggregate into droplets with free interfaces that separate the high density nematic from the low density isotropic phases, depicted in Fig. 1(iv). This is conveniently described by a nondimensional “phase field” corresponding to filament density, $\psi$, where $\psi > 0$ ($\psi < 0$) describes the interior (exterior) of the droplet. Filaments in the high density droplet interior align in orientation, described by the 2D nematic order parameter, $Q_{ij}$. These two ingredients result in nematic droplets with tactoid shape at equilibrium [36]. The active motion of filaments and polarity sorting induced by motors results in a net polarization within the droplet. Unlike in the 1D model, polar order in 2D is non-conserved and can be induced by motor-driven torques or relaxed by rotational diffusion. Generalizing Eq. (1) - (3) to 2D, and observing usual principles of conservation and symmetry [5], we obtain the dynamical equations,

$$ \partial_t m = D_m \nabla^2 m + v_0 \nabla \cdot (mp) + k_{\text{on}} \psi \theta(\psi) - k_{\text{off}} m, \quad (4) $$

$$ \partial_t \psi = \Gamma_{\psi} \nabla^2 \frac{\delta F}{\delta \psi} - \zeta \nabla \cdot (mp), \quad (5) $$

$$ \partial_t p = -\zeta_0 \nabla \psi - \zeta_p \nabla m - \Gamma_p \frac{\delta F}{\delta p}, \quad (6) $$

$$ \partial_t Q_{ij} = -\Gamma_Q \frac{\delta F}{\delta Q_{ij}}. \quad (7) $$

Eq. (4)-(5) describe the conservation of motors and filaments, respectively, and include both active and passive fluxes on the right side. The motor flux includes diffusion in the first term, and active motion of motors with a propulsion velocity $v_0$ in the second. Additionally, we include binding-unbinding kinetics for the motors (third and fourth term in Eq. (4)), where motors bind with rate $k_{\text{on}}$ wherever filaments exist (expressed by $\theta(\psi)$, the Heaviside step function) from a “reservoir” of free motors in solution, and unbind with rate $k_{\text{off}}$, but are not restricted to diffuse within the droplet. Eq. (5) includes a flux created by a free energy of inter-filament interactions (full form given below) as well as active flux induced
by motors advecting filaments in their polarity direction. Eq. (6) describes the induction of polarization by torques caused by gradients in motor and filament density, and its relaxation through rotational diffusion. Seen in the framework of the Toner-Tu hydrodynamic theory that describes the flocking of polar active particles [37], where p has the status of both an orientational order parameter as well as a fluid velocity, the $-\zeta_0 \nabla \psi$ gives the gradient of pressure, and $-\zeta_p \nabla m$ is the gradient of an active stress created by motors [38]. Terms similar to this latter also arise in the theory of chemotactic colloids [39] and equilibrium polar liquid crystals [40], in both of which cases a chemical concentration can guide polarization. Microscopically, the $\zeta_p$ term captures the preferred orientation of filaments towards regions with higher motor density by the binding and crosslinking by motors, while $\zeta$ originates from the sliding forces exerted by motors. Note that the $\zeta$ specifies both polarization and active flux in the one dimensional model, whereas the $\zeta_p$ term arises in the two dimensional model because of the additional orientational degree of freedom. The parameters $v_0$, $\zeta$ and $\zeta_p$ then all depend on motor activity but can in principle be varied independently by tuning motor properties such as size, shape, processivity and crosslinking. We assume in Eq. (7) that the nematic order is strong and arises from equilibrium forces. The timescales for relaxation towards equilibrium are specified by the “mobility” coefficients, $\Gamma_\psi$, $\Gamma_p$ and $\Gamma_Q$.

The equilibrium dynamics assume a coupled phase transition in $\psi$ and nematic order $Q$, and a relaxation of $p$ included in the total free energy,

$$F = F_{\psi} + F_{p} + F_{Q} + F_{\text{int}}, \quad (8a)$$

$$F_{\psi} = -\frac{\nu_2}{2}\psi^2 + \frac{\nu_4}{4}\psi^4 + \frac{B}{2}(\nabla\psi)^2, \quad (8b)$$

$$F_p = \frac{\alpha_p}{2} |\mathbf{p}|^2 + \frac{\beta_p}{2} |\mathbf{p}|^4 + \frac{\kappa_p}{2} (\nabla \cdot \mathbf{p})^2, \quad (8c)$$

$$F_Q = -\frac{\alpha_Q}{4} Q_{ij}^2 + \frac{\beta_Q}{16} Q_{ij}^2 Q_{kh}^2 + \frac{\kappa_Q}{2} (\partial_k Q_{ik})^2, \quad (8d)$$

$$F_{\text{int}} = \frac{C}{4} Q_{ij}^2 \psi(\psi - 1) + A \partial_i \psi \partial_j Q_{ij} + \Omega Q_{ij} p_i p_j. \quad (8e)$$

The density free energy Eq. (8b) models phase separating droplets according to standard Cahn-Hilliard dynamics. Ignoring corrections for curved interfaces [41], we use droplet surface tension (“line tension” in 2D) as $\gamma = \frac{2\sqrt{2B\nu_2}}{3\nu_4}$ and its interfacial width, $\xi = \sqrt{2B/\nu_2}$ [42]. The free energy for the polarization, Eq. (8c), includes two relaxation terms, $\alpha_p$, $\beta_p > 0$ (corresponding to lack of spontaneous polar order) and an elastic term, $\kappa_p$. Eq. (8d) is the Landau-de Gennes free energy for the nematic order, $Q_{ij}$, with elasticity $\kappa_Q$. All equilibrium couplings between fields are written in Eq. (8e), where the first term controls the density-driven isotropic-nematic transition and induces nematic order within the droplet. The second term is a “weak anchoring” that aligns the nematic parallel to the droplet interface for $A > 0$. The third term aligns the polarization with the nematic or-

<figure><img src="image_3.png" /><figcaption>FIG. 2. (i) and (ii) A representative numerical steady state solution of our active droplet model. (i) Red gradient indicates the density $\psi$ and the white lines indicate the nematic director field. (ii) Blue gradient, black arrows and the red line indicate the motor concentration, the polarization direction and the droplet boundary respectively. (iii): Time series for different droplet types observed in our simulations: (iii a) single tactoid with centering aster; (iii b) two connected tactoids with centering aster; (iii c) three connected tactoids with two asters; (iii d) fully divided droplet. The number in each snapshot shows the simulation time.</figcaption></figure>

der. Overall, this model recapitulates the key elements of nematic phase separation of the filaments and their coupling with motor activity.

**Results.** We now employ numerical simulations to explore the consequences of motor activity on the dynamics and morphology of phase-separating nematic droplets. Starting from an initially circular droplet of radius $R_0$, we integrate Eqs. 4-7 until a non-equilibrium steady-state is reached (see SI III for simulation details and parameters). Figure 2 presents typical simulation results. We obtain an elongated droplet with the nematic aligned along its long axis, as shown by the density and nematic director plots for a typical case in Fig. 2(i). The motor density accumulates at the core of the aster it induces with an outward polarization as shown in Fig. 2(ii). We first explore the interplay between motor-generated active forces, which tend to distort the equilibrium structure, and surface tension which resists such deformation. To this aim, we perform simulations with varying mo-
<figure><img src="image_4.png" /><figcaption>FIG. 3. (i) Nonequilibrium phase diagram from our simulations for varying surface tension $\gamma$ and $\zeta_p$ while keeping $\zeta = 1.5$ fixed. The solid line is a fit to our theoretical analysis. (ii)-(v) Different steady state shapes observed in our simulations and experiments (scale bars are 5µm). For the simulation results (ii)-(v) (top) colorcode is the motor concentration $c$ as, black arrows are the polarization, $p$, and the red contour shows the shape of the droplet. We classify the droplet shapes as follows: (ii) single tactoid with a centering aster; (iii) two tactoids with one aster; (iv) three connected tactoids with two asters; (v) fully divided droplet. (vi) Critical motor concentration for the transition from one to two asters as a function of initial droplet size. (vii) Nonequilibrium phase diagram showing the different steady states observed in our simulations for varying activity $\zeta$ and surface tension $\gamma$ while keeping $\zeta_p = 5.0$ fixed. The solid line is a fit to the prediction from our theoretical analysis.</figcaption></figure>

tor induced polarization, $\zeta_p$, relative to surface tension, $\gamma$. The resulting time sequences of observed droplet shapes, starting from an unpolarized droplet, are shown in Fig. 2(iii). When $\zeta_p$ is low compared to surface tension, Fig. 2(iii a), motors localize towards the center to form an aster, which only slightly polarizes and deforms the droplet. At intermediate $\zeta_p$ (iii b), motors localize more strongly, resulting in a stronger elongation and the appearance of a constriction of the droplet around the midplane between strongly polarized lobes. Increasing the $\zeta_p$ further can lead to two distinct scenarios: In (iii c), the aster divides and a third, central lobe with strong polarity gradient emerges between two constrictions. In (iii d), finally, $\zeta_p$ is sufficiently strong to induce full division of the initial droplet into two polarized daughter droplets. The short-time dynamics of this model thus reproduces the motor centering, aster formation and polarity sorting features of the 1D model (Fig. 1), while at longer times a complex diversity of morphologies emerges from the interplay between surface tension and motor activity.

To rationalize this rich phenomenology, we build a morphological phase diagram in Fig. 3(i) in $\zeta_p-\gamma$ while keeping other parameters constant. Each distinct steady

state structure is shown in Fig. 3(ii)-(v) (top). Interestingly, each of these morphologies corresponds to experimentally observed shapes, as shown in Fig. 3(ii)-(v) (bottom). Confirming the qualitative findings described in Fig. 2(iii), at medium to high surface tension and low $\zeta_p$, we find motors form asters in the midplane of the undeformed droplet (blue stars), whereas higher $\zeta_p$ increases the influence of the centered aster on droplet shape, pinching it into two lobes (orange circles). This is consistent with the experimental observation that motors always localize to droplet midplane, but only deform the droplet when there are more active motors [31]. Qualitatively, the motors at the aster core splay the filaments, which are anchored to the interface, and which therefore results in an inward pinching of the interface. The motor induced polarization (second term in Eq. 6) is derivable from a functional, $-(\zeta_p/\Gamma_p) \int d^2x m \partial_i p_i$, corresponding to an effective free energy for motor-induced spontaneous splay. This corresponds to a negative contribution to surface energy [43] (see SI IV for the derivation), that drives the droplet shape instability when $\zeta_p m_0 / \Gamma_p > \gamma$. We thus expect the transition between undeformed and pinched droplets to occur across a dividing line, $\zeta_p \sim \gamma$, which is confirmed by our phase diagram data in Fig. 3(i).
By further increasing $\zeta_p$ and lowering surface tension, we find droplets with two asters which pinch them into three-lobed structures (Fig. 3(iv), green upward triangles). Importantly, we also find experimental realizations of such multiply pinched, equi-lobed structures, shown in Fig. 3(iv) bottom panel. Such linear “strings of tactoids” connected by multiple motor clusters have not been previously reported and are evocative of fibers with periodic contractile units such as in muscle or anomalous, multi-polar biological spindles [44].

Multi-aster states such as aster lattices are in fact a generic feature of bulk active polar fluids [26, 30], but here, we report and analyze their occurrence within droplets. To explore the aster-forming instability arising from the feedback between motor flux and motor induced polarization in Eq. 4 and Eq. 6, we perform a linear stability analysis for an incompressible, polar fluid featuring only $p$ and $m$ (see SI V). The analysis yields a characteristic spacing between asters in bulk, $l_c$, such that $l_c \sim 1/\sqrt{m_0\zeta_p}$, that decreases with motor density and motor induced polarization. Based on this, we expect multiple asters can be accommodated within the droplet when its major axis, $l_d$, is long enough compared to $l_c$. We confirm this numerically in (Fig. 3(vi)) by showing that the critical motor concentration for transition from one to two asters for varying droplet size (but all other parameters constant), does indeed decrease with initial droplet size, with an expected scaling of $m_0 \sim 1/l_d^2$. Since the aspect ratio of a droplet at equilibrium increases with decreasing surface tension, the droplets get longer and can accommodate multiple asters when $l_d > l_c \sim \zeta_p^{-1/2}$ is satisfied. This explains why lower surface tension favors the two-aster state (green triangles) in Fig. 3(i). At even higher $\zeta_p$ and lower $\gamma$, we find a region of the phase diagram in Fig. 3(i) where the droplet is fully divided (purple diamonds).

To see how droplet division can be enhanced, we now explore an alternative shape instability in our model Eq. 5, whereby an active density flux pushes filaments away from the aster. By varying the strength of this active flux, $\zeta$, and the surface tension $\gamma$, while keeping $\zeta_p(= 5.0)$ fixed, we obtain the phase diagram in Fig. 3(vii), which shows undeformed single tactoids (blue stars), pinched droplets (orange circles) and fully divided droplets (purple diamonds). Here, we find that droplet division occurs over a wide range in parameter space, showing that $\zeta$ is a good control parameter to trigger droplet division, while $\zeta_p$ can be used to obtain multi-asters states.

To analyze how the active flux parameter, $\zeta$, for full division scales with $\gamma$, we note that if the motor and polarization are “fast variables” that relax quickly to their steady state aster solution (shown in the SI VI), the active density flux term in Eq. (5), $-\zeta\nabla \cdot (mp)$ can be

obtained from the variation of an effective “free energy” functional, $-\zeta \int d^2x \omega\psi$. Here, $\omega$ is a scalar potential that gives the steady state polarization corresponding to an aster through $mp = -\nabla\omega$. It can be shown that this effective free energy term makes a negative contribution to the droplet surface energy (detailed in the SI VI) that scales with $\zeta$ and which can then destabilize the droplet when sufficiently strong compared to surface tension. We thus expect the transition to fully divided droplets to scale as $\zeta \sim \gamma$, which is confirmed by our simulations in Fig. 3(vii).

*Conclusion.* We build and explore a minimal theoretical model for motor-filament droplets which captures four different shapes and show them experimentally in the actomyosin droplet system. We predict a phase diagram of expected droplet shapes based on motor activity and filament interactions, which can be tested in future experiments by systematically varying these properties. Unlike other theoretically proposed active droplet division mechanisms, the droplet shape changes we predict and show in experiment, do not require fluxes of chemical reactants [13] or large-scale fluid flows and defect dynamics [12], that arise in active nematic morphodynamics [45]. Instead, the emergence of local polar order within a nematic and the consequent spatial localization of activity, is crucial in our model, which was not explored in a previous equilibrium nematic model for this system [31]. In addition to *in vitro* cytoskeletal spindles, our results also have implications for the physical principles behind cell division and the organization of the biological spindle, where motor-induced active organization of filaments [9], shape change [46], co-existence of polar order with nematic order [47] and even formation of multi-lobed structures [48] can occur. We expect our work to inform strategies to use localization of active agents to achieve self-actuated shape morphing in materials [49] and synthetic cells.

## ACKNOWLEDGMENTS

FJS and KD gratefully acknowledge computing time on the Multi-Environment Computer for Exploration and Discovery (MERCED) cluster at the University of California, Merced, which was funded by National Science Foundation Grant No. ACI-1429783 and support from the National Science Foundation: NSF-CREST: Center for Cellular and Biomolecular Machines (CCBM) at the University of California, Merced: NSF-HRD-1547848.. PR received funding from the “Investissements d’Avenir” French Government program managed by the French National Research Agency (ANR-16-CONV-0001) and from Excellence Initiative of Aix-Marseille University - A*MIDEX.

<footnote>[1] M. Murrell, P. W. Oakes, M. Lenz, and M. L. Gardel, Forcing cells into shape: the mechanics of actomyosin</footnote>
486 EP (2015), review Article.

[2] V. Noireaux, Y. T. Maeda, and A. Libchaber, Development of an artificial cell, from self-organization to computation and self-reproduction, Proceedings of the National Academy of Sciences 108, 3473 (2011), https://www.pnas.org/content/108/9/3473.full.pdf.

[3] M. Soares e Silva, M. Depken, B. Stuhrmann, M. Korsten, F. C. MacKintosh, and G. H. Koenderink, Active multistage coarsening of actin networks driven by myosin motors, Proceedings of the National Academy of Sciences of the United States of America 108, 9408 (2011), 21593409[pmid].

[4] T. Sanchez, D. T. Chen, S. J. DeCamp, M. Heymann, and Z. Dogic, Spontaneous motion in hierarchically assembled active matter, Nature 491, 431 (2012).

[5] M. C. Marchetti, J. F. Joanny, S. Ramaswamy, T. B. Liverpool, J. Prost, M. Rao, and R. A. Simha, Hydrodynamics of soft active matter, Rev. Mod. Phys. 85, 1143 (2013).

[6] A. A. Hyman, C. A. Weber, and F. Jülicher, Liquid-liquid phase separation in biology, Annual Review of Cell and Developmental Biology 30, 39 (2014), pMID: 25288112, https://doi.org/10.1146/annurev-cellbio-100913-013325.

[7] P.-G. De Gennes and J. Prost, The Physics of Liquid Crystals, Vol. 83 (Oxford University Press, 1995).

[8] K. L. Weirich, S. Banerjee, K. Dasbiswas, T. A. Witten, S. Vaikuntanathan, and M. L. Gardel, Liquid behavior of cross-linked actin bundles, Proceedings of the National Academy of Sciences 114, 2131 (2017).

[9] J. Brugués and D. Needleman, Physical basis of spindle self-organization, Proceedings of the National Academy of Sciences 111, 18496 (2014), https://www.pnas.org/content/111/52/18496.full.pdf.

[10] A. V. Kaznacheev, M. M. Bogdanov, and S. A. Taraskin, The nature of prolate shape of tactoids in lyotropic inorganic liquid crystals, Journal of Experimental and Theoretical Physics 95, 57 (2002).

[11] P. Prinsen and P. van der Schoot, Shape and director-field transformation of tactoids, Phys. Rev. E 68, 021701 (2003).

[12] L. Giomi and A. DeSimone, Spontaneous division and motility in active nematic droplets, Phys. Rev. Lett. 112, 147802 (2014).

[13] D. Zwicker, R. Seyboldt, C. A. Weber, A. A. Hyman, and F. Jülicher, Growth and division of active droplets provides a model for protocells, Nature Physics 13, 408 (2017).

[14] H. R. Vutukuri, M. Hoore, C. Abaurrea-Velasco, L. van Buren, A. Dutto, T. Auth, D. A. Fedosov, G. Gompper, and J. Vermant, Active particles induce large shape deformations in giant lipid vesicles, Nature 586, 52 (2020).

[15] C. Wang, Y.-k. Guo, W.-d. Tian, and K. Chen, Shape transformation and manipulation of a vesicle by active particles, The Journal of chemical physics 150, 044907 (2019).

[16] Y. Li and P. R. ten Wolde, Shape transformations of vesicles induced by swim pressure, Physical review letters 123, 148003 (2019).

[17] S. C. Takatori and A. Sahu, Active contact forces drive nonequilibrium fluctuations in membrane vesicles, Physical Review Letters 124, 158102 (2020).

[18] R. Singh and M. Cates, Hydrodynamically interrupted droplet growth in scalar active matter, Physical review letters 123, 148005 (2019).

[19] T. Ruiz-Herrero, T. G. Fai, and L. Mahadevan, Dynamics of growth and form in prebiotic vesicles, Physical review letters 123, 038102 (2019).

[20] E. Tjhung, D. Marenduzzo, and M. E. Cates, Spontaneous symmetry breaking in active droplets provides a generic route to motility, Proceedings of the National Academy of Sciences 109, 12381 (2012), https://www.pnas.org/content/109/31/12381.full.pdf.

[21] F. Ziebert, S. Swaminathan, and I. S. Aranson, Model for self-polarization and motility of keratocyte fragments, Journal of The Royal Society Interface 9, 1084 (2012).

[22] F. C. Keber, E. Loiseau, T. Sanchez, S. J. DeCamp, L. Giomi, M. J. Bowick, M. C. Marchetti, Z. Dogic, and A. R. Bausch, Topology and dynamics of active nematic vesicles, Science 345, 1135 (2014).

[23] T. Surrey, F. Nédélec, S. Leibler, and E. Karsenti, Physical properties determining self-organization of motors and microtubules, Science 292, 1167 (2001).

[24] H. Youn Lee and M. Kardar, Macroscopic equations for pattern formation in mixtures of microtubules and molecular motors, Phys. Rev. E 64, 056113 (2001).

[25] K. Kruse, J. F. Joanny, F. Jülicher, J. Prost, and K. Sekimoto, Asters, vortices, and rotating spirals in active gels of polar filaments, Phys. Rev. Lett. 92, 078101 (2004).

[26] S. Sankararaman, G. I. Menon, and P. B. Sunil Kumar, Self-organized pattern formation in motor-microtubule mixtures, Phys. Rev. E 70, 031905 (2004).

[27] I. S. Aranson and L. S. Tsimring, Theory of self-assembly of microtubules and motors, Phys. Rev. E 74, 031915 (2006).

[28] A. Ahmadi, T. B. Liverpool, and M. C. Marchetti, Nematic and polar order in active filament solutions, Phys. Rev. E 72, 060901 (2005).

[29] F. Ziebert and W. Zimmermann, Nonlinear competition between asters and stripes in filament-motor systems, The European Physical Journal E 18, 41 (2005).

[30] K. Gowrishankar and M. Rao, Nonequilibrium phase transitions, fluctuations and correlations in an active contractile polar fluid, Soft matter 12, 2040 (2016).

[31] K. L. Weirich, K. Dasbiswas, T. A. Witten, S. Vaikuntanathan, and M. L. Gardel, Self-organizing motors divide active liquid droplets, Proceedings of the National Academy of Sciences 116, 11125 (2019).

[32] K. Kruse and F. Jülicher, Self-organization and mechanical properties of active filament bundles, Phys. Rev. E 67, 051913 (2003).

[33] M. Stachowiak, P. McCall, T. Thoresen, H. Balcioglu, L. Kasiewicz, M. Gardel, and B. O'Shaughnessy, Self-organization of myosin ii in reconstituted actomyosin bundles, Biophysical Journal 103, 1265 (2012).

[34] N. Kumar, R. Zhang, J. J. de Pablo, and M. L. Gardel, Tunable structure and dynamics of active liquid crystals, Science Advances 4, 10.1126/sciadv.aat7779 (2018).

[35] S. Fürthauer, B. Lemma, P. J. Foster, S. C. EmsMcClung, C.-H. Yu, C. E. Walczak, Z. Dogic, D. J. Needleman, and M. J. Shelley, Self-straining of actively crosslinked microtubule networks, Nature Physics 15, 1295 (2019).

[36] N. B. Ludwig, K. L. Weirich, E. Alster, T. A. Witten, M. L. Gardel, K. Dasbiswas, and S. Vaikuntanathan, Nucleation and shape dynamics of model nematic tactoids around adhesive colloids, The Journal of Chemical Physics 152, 084901 (2020), https://doi.org/10.1063/1.5141997.
[37] J. Toner and Y. Tu, Long-range order in a two-dimensional dynamical XY model: How birds fly together, Phys. Rev. Lett. 75, 4326 (1995).

[38] K. Husain and M. Rao, Emergent structures in an active polar fluid: Dynamics of shape, scattering, and merger, Phys. Rev. Lett. 118, 078104 (2017).

[39] B. Liebchen and H. Löwen, Synthetic chemotaxis and collective behavior in active matter, Accounts of chemical research 51, 2982 (2018).

[40] W. Kung, M. Cristina Marchetti, and K. Saunders, Hydrodynamics of polar liquid crystals, Phys. Rev. E 73, 031708 (2006).

[41] J. C. Everts, M. T. J. J. M. Punter, S. Samin, P. van der Schoot, and R. van Roij, A landau-de gennes theory for hard colloidal rods: Defects and tactoids, The Journal of Chemical Physics 144, 194901 (2016), https://doi.org/10.1063/1.4948785.

[42] S. Safran, Statistical thermodynamics of surfaces, interfaces, and membranes (CRC Press, 2018).

[43] D. Pettey and T. C. Lubensky, Stability of texture and shape of circular domains of langmuir monolayers, Phys. Rev. E 59, 1834 (1999).

[44] N. J. Ganem, S. A. Godinho, and D. Pellman, A mechanism linking extra centrosomes to chromosomal instability, Nature 460, 278 (2009).

[45] L. Metselaar, J. M. Yeomans, and A. Doostmohammadi, Topology and morphology of self-deforming active shells, Phys. Rev. Lett. 123, 208001 (2019).

[46] D. Oriola, F. Jülicher, and J. Brugués, Active forces shape the metaphase spindle through a mechanical instability, Proceedings of the National Academy of Sciences 117, 16154 (2020), https://www.pnas.org/content/117/28/16154.full.pdf.

[47] J. Roostalu, J. Rickman, C. Thomas, F. Nédélec, and T. Surrey, Determinants of polar versus nematic organization in networks of dynamic microtubules and mitotic motors, Cell 175, 796 (2018).

[48] M. Leoni, O. V. Manyuhina, M. J. Bowick, and M. C. Marchetti, Defect driven shapes in nematic droplets: analogies with cell division, Soft Matter 13, 1257 (2017).

[49] D. P. Holmes, Elasticity and stability of shape-shifting structures, Current opinion in colloid & interface science 40, 118 (2019), http://doi.org/10.1016/j.cocis.2019.02.008.

[50] J. Liouville, Sur l'équation aux différences partielles $\frac{d^2 \log \lambda}{dudv} \pm \frac{\lambda}{2a^2} = 0$ $\frac{d^2 \log \lambda}{dudv} \pm \frac{\lambda}{2a^2} = 0$, Journal de mathématiques pures et appliquées, 71 (1853).

[51] G. Bratu, Sur les équations intégrales non linéaires, Bulletin de la Société Mathématique de France 42, 113-142 (1914).

[52] I. M. Gelfand, Some problems in the theory of quasilinear equations, Amer. Math. Soc. Transl 29.2, 295-381 (1963).

[53] G. W. Walker, Some problems illustrating the forms of nebulae, Proceedings of the Royal Society of London. Series A, Containing Papers of a Mathematical and Physical Character 91.631, 410 (1915).
<figure><img src="image_5.png" /><figcaption>FIG. S1. Division of a one dimensional droplet. Each plot shows a time snapshot with the time shown in the upper left corner. Blue solid lines show the motor concentration, orange dashed lines show the polarization and green dashed dotted lines show the density.</figcaption></figure>

# Supplementary Information

## I. DERIVATION AND SOLUTION OF ONE DIMENSIONAL MODEL

The 1D model for the dynamics of filament density, filament polarization and motor density given in the main text follow from the conservation equations:

$$ \partial_t n_+ = D \partial_x^2 n_+ - \zeta \partial_x (n_+ m), \qquad (S1) $$

$$ \partial_t n_- = D \partial_x^2 n_- + \zeta \partial_x (n_- m), \qquad (S2) $$

$$ \partial_t m = D_m \partial_x^2 m + v_0 \partial_x m (n_+ - n_-), \qquad (S3) $$

Our one dimensional model (Eqs. (1)-(3)) is simplified by assuming incompressibility, i.e. the density is $\rho = \rho_0$. In the steady state Eqs. (1)-(3) then reduce to

$$ \partial_t p = D \partial_x^2 p - \zeta \rho_0 \partial_x m, \qquad (S4) $$

$$ \partial_t m = D_m \partial_x^2 m + v_0 \partial_x (m p). \qquad (S5) $$

The Eqs. (S4)-(S5) are solved with “natural boundary”conditions that the motor density and polarization as well as their corresponding fluxes decay to zero at $x = \pm \infty$ to yield:

$$ p(x) = \tanh(x/\xi), \qquad (S6) $$

$$ m(x) = m_0 \mathrm{sech}^2(x/\xi), \qquad (S7) $$

where $\xi = \frac{2D_m}{v_0}$ is the typical aster length scale and $m_0 = \frac{D\xi}{\zeta\rho_0}$ is the aster strength. Figure 1(ii-iii) shows the solutions Eqs. (S6)-(S7), which show the one dimensional equivalent of an out pointing aster. The motors gather at the center of the aster as expected. Note that in 1D, the flux at steady state is a constant throughout the system, and that natural boundary conditions require this flux to vanish. In other words the diffusive flux of the motors and filaments is balanced by the corresponding active flux that advects them, at steady state. This competition sets the width of the distribution of motors at the “aster”: $\xi \sim \frac{D_m}{v_0}$. In a confined droplet geometry, and at higher dimensions, solutions with nonzero fluxes are possible and seen in our numeric solutions.

## II. 1D NUMERICS

To test the minimal 1D model in a finite domain, we numerically integrated the following equations

$$ \partial_t \rho = D \partial_x^2 \rho - \zeta \rho_0 \partial_x (p m) + \Gamma_{\rho} \partial_x^2 \frac{\delta F_{1d}}{\delta \rho}, \qquad (S8) $$

$$ \partial_t p = D \partial_x^2 p - \zeta / \rho_0 \partial_x (\rho m) - \Gamma_p \frac{\delta F_{1d}}{\delta p}, \qquad (S9) $$

$$ \partial_t m = D_m \partial_x^2 m + v_0 \partial_x (m p). \qquad (S10) $$
<figure><img src="image_6.png" /><figcaption>FIG. S2. Single droplet with centering aster in one dimensional model. (a) Motor concentration. (b) Polarization o filaments. (c) Density of filaments.</figcaption></figure>

Here, we included the following free energy

$$F_{1d} = -\frac{\nu_2}{2}\rho^2 + \frac{\nu_4}{4}\rho^4 + \frac{B}{2}(\partial_x\rho)^2 \quad (S11)$$

$$+ \alpha_p p^2/2 + \beta_p p^4/2, \quad (S12)$$

which has a phase field for the density $\rho$ such that a one dimensional droplet is formed and $F_{1d}$ has decaying terms for the polarization $p$.

We integrate Eq. (S8)-(S10) numerically using a finite difference method with a fourth order Runge-Kutta time discretization and periodic boundary conditions. The simulation parameters can be found in table I.

Figure S2 show the motor concentration (Fig. S2(a)), polarization (Fig. S2(b)) and density (Fig. S2(c)) for a simulation at low activity $\zeta = 1.5$. We find the emergence of a single droplet with a centering aster.

Going to higher activity $\zeta = 3.5$, we find that the droplet can divide by means of the active motion of filaments. Figure S1 shows a time series of the division of a droplet. Here, the blue solid line show the motor concentration, the dashed orange line the polarization, and the dashed dotted green line the density.

### III. SIMULATION METHOD

We solve Eq. (4)-(7) numerically using a finite difference scheme with a fourth order Runge-Kutta time integration on a 74 × 74 grid. We use a timestep of $\delta t = 0.01$ and the various parameters are given in Table III, while we vary surface tension, $\gamma$, and the strength of motor-induced density flux and molecular torque: $\zeta$ and $\zeta_p$. We initialize the system with a circular droplet of radius $R_0$.

The parameters used in the simulations in the main text are given in Table II and Table III.
<table><caption>TABLE I. Parameter values used in one dimensional model.</caption><thead><tr><th>Parameter</th><th>Value</th></tr></thead><tbody><tr><td>Δx</td><td>1.0</td></tr><tr><td>grid</td><td>128</td></tr><tr><td>δt</td><td>0.001</td></tr><tr><td>D<sub>m</sub></td><td>2</td></tr><tr><td>v<sub>0</sub></td><td>2</td></tr><tr><td>ζ</td><td>10.0</td></tr><tr><td>κ<sub>p</sub></td><td>4.0</td></tr><tr><td>α<sub>p</sub></td><td>0.1</td></tr><tr><td>β<sub>p</sub></td><td>1.0</td></tr><tr><td>Γ<sub>p</sub></td><td>1.0</td></tr><tr><td>ν<sub>2</sub></td><td>0.01</td></tr><tr><td>ν<sub>4</sub></td><td>1.4</td></tr><tr><td>B</td><td>0.04</td></tr><tr><td>Γ<sub>ψ</sub></td><td>1.0</td></tr><tr><td>R<sub>0</sub></td><td>35</td></tr><tr><td>c<sub>0</sub></td><td>0.1</td></tr></tbody></table>

<table><caption>TABLE II. Parameter values used for the example droplet shapes in Fig. 2</caption><thead><tr><th>Parameter</th><th>Fig. 2iiia</th><th>Fig.2(c)</th><th>Fig.2(h)</th><th>Fig.2(j)</th></tr></thead><tbody><tr><td>Δx Δy</td><td>1.0 1.0</td><td>1.0 1.0</td><td>1.0 1.0</td><td>1.0 1.0</td></tr><tr><td>grid</td><td>100 × 100</td><td>100 × 100</td><td>100 × 100</td><td>100 × 100</td></tr><tr><td>D<sub>m</sub></td><td>0.4</td><td>0.4</td><td>0.4</td><td>0.4</td></tr><tr><td>v<sub>0</sub></td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td></tr><tr><td>k<sub>on</sub></td><td>0.00005</td><td>0.00015</td><td>0.00015</td><td>0.00015</td></tr><tr><td>k<sub>off</sub></td><td>0.0005</td><td>0.001</td><td>0.001</td><td>0.001</td></tr><tr><td>ν<sub>2</sub></td><td>0.4</td><td>0.4</td><td>0.4</td><td>0.4</td></tr><tr><td>ν<sub>4</sub></td><td>2.0</td><td>2.0</td><td>2.0</td><td>2.0</td></tr><tr><td>Γ<sub>ψ</sub></td><td>1.0</td><td>1.0</td><td>1.0</td><td>1.0</td></tr><tr><td>R<sub>0</sub></td><td>23.0</td><td>23.0</td><td>23.0</td><td>23.0</td></tr><tr><td>κ<sub>p</sub></td><td>0.1</td><td>0.1</td><td>0.1</td><td>0.1</td></tr><tr><td>α<sub>p</sub></td><td>0.1</td><td>0.1</td><td>0.1</td><td>0.1</td></tr><tr><td>β<sub>p</sub></td><td>10.0</td><td>10.0</td><td>10.0</td><td>10.0</td></tr><tr><td>Γ<sub>p</sub></td><td>1.0</td><td>1.0</td><td>1.0</td><td>1.0</td></tr><tr><td>ζ<sub>0</sub></td><td>0.01</td><td>0.01</td><td>0.01</td><td>0.01</td></tr><tr><td>κ<sub>Q</sub></td><td>0.55</td><td>0.4</td><td>0.4</td><td>0.4</td></tr><tr><td>α<sub>Q</sub></td><td>0.01</td><td>0.01</td><td>0.01</td><td>0.01</td></tr><tr><td>β<sub>Q</sub></td><td>2.0</td><td>2.0</td><td>2.0</td><td>2.0</td></tr><tr><td>Γ<sub>Q</sub></td><td>1.0</td><td>1.0</td><td>1.0</td><td>1.0</td></tr><tr><td>C</td><td>1.0</td><td>1.0</td><td>1.0</td><td>1.0</td></tr><tr><td>A</td><td>0.6</td><td>0.8</td><td>0.8</td><td>0.8</td></tr><tr><td>Ω</td><td>0.3</td><td>0.3</td><td>0.3</td><td>0.3</td></tr><tr><td>B</td><td>1.8</td><td>2.8</td><td>2.8</td><td>2.8</td></tr><tr><td>ζ</td><td>4.0</td><td>1.0</td><td>0.8</td><td>1.0</td></tr><tr><td>ζ<sub>p</sub></td><td>20.0</td><td>40.0</td><td>5.0</td><td>20.0</td></tr></tbody></table>

## IV. ASTER CONTRIBUTION TO SURFACE ENERGY

The formation of an aster is actively induced by the motor gradient term in Eq. 6, which can be shown to arise from an effective free energy contribution, $ζ_p/Γ_p \int d^2x \mathbf{p} \cdot \nabla m$, to the polarization equation. Integrating by parts and setting the surface term to zero, this term can be seen to be equivalent to a spontaneous splay energy, $-ζ_p/Γ_p \int d^2x m \nabla \cdot \mathbf{p}$.

For a uniform motor concentration, $m(x) = m_0$, this spontaneous splay term can be transformed into a surface term using the Green's theorem (which is also the 2D version of the divergence theorem) to give, $-ζ_p m_0/Γ_p \int dl \hat{\mathbf{n}} \cdot \mathbf{p}$, where the line element $dl$ is along the boundary of the droplet and $\hat{\mathbf{n}}$ is the unit normal to the droplet boundary. Since the $\mathbf{p}$ points outwards in the asters, we get $\hat{\mathbf{n}} \cdot \mathbf{p} > 0$, which therefore contributes a negative line tension to the droplet interfacial energy. Comparing with the droplet line tension energy, $\gamma \int dl$, and assuming that polarization has a value $p_0$ at the droplet boundary, we arrive at the condition, $ζ_p p_0 m_0 > Γ_p \gamma$ for destabilization of droplet shape by
<table><caption>TABLE III. Parameter values used in the simulations to compute the phase diagrams (Fig.3(i) and Fig.3(vii)) and to compute the critical motor concentration for three aster splitting as function of initial droplet size (Fig.3(vi)).</caption><thead><tr><th colspan="2">Fig.3(i),(vii))</th><th colspan="2">Fig.3(vi))</th></tr><tr><th>Parameter</th><th>Value</th><th>Parameter</th><th>Value</th></tr></thead><tbody><tr><td>Δx Δy</td><td>1.0 1.0</td><td>Δx Δy</td><td>1.0 1.0</td></tr><tr><td>grid</td><td>74 × 74</td><td>grid</td><td>74 × 74</td></tr><tr><td>D<sub>m</sub></td><td>0.4</td><td>D<sub>m</sub></td><td>0.4</td></tr><tr><td>v<sub>0</sub></td><td>0.15</td><td>v<sub>0</sub></td><td>0.15</td></tr><tr><td>k<sub>on</sub></td><td>0.00015</td><td>k<sub>off</sub></td><td>0.001</td></tr><tr><td>k<sub>off</sub></td><td>0.001</td><td>ν<sub>2</sub></td><td>0.4</td></tr><tr><td>ν<sub>2</sub></td><td>0.4</td><td>ν<sub>4</sub></td><td>2.0</td></tr><tr><td>ν<sub>4</sub></td><td>2.0</td><td>Γ<sub>ψ</sub></td><td>1.0</td></tr><tr><td>Γ<sub>ψ</sub></td><td>1.0</td><td>κ<sub>p</sub></td><td>0.1</td></tr><tr><td>R<sub>0</sub></td><td>17.0</td><td>α<sub>p</sub></td><td>0.1</td></tr><tr><td>κ<sub>p</sub></td><td>0.1</td><td>β<sub>p</sub></td><td>10.0</td></tr><tr><td>α<sub>p</sub></td><td>0.1</td><td>Γ<sub>p</sub></td><td>1.0</td></tr><tr><td>β<sub>p</sub></td><td>10.0</td><td>ζ<sub>0</sub></td><td>0.01</td></tr><tr><td>Γ<sub>p</sub></td><td>1.0</td><td>κ<sub>Q</sub></td><td>0.3</td></tr><tr><td>ζ<sub>0</sub></td><td>0.01</td><td>α<sub>Q</sub></td><td>0.01</td></tr><tr><td>κ<sub>Q</sub></td><td>0.3</td><td>β<sub>Q</sub></td><td>2.0</td></tr><tr><td>α<sub>Q</sub></td><td>0.01</td><td>Γ<sub>Q</sub></td><td>1.0</td></tr><tr><td>β<sub>Q</sub></td><td>2.0</td><td>C</td><td>1.0</td></tr><tr><td>Γ<sub>Q</sub></td><td>1.0</td><td>A</td><td>0.5</td></tr><tr><td>C</td><td>1.0</td><td>Ω</td><td>0.3</td></tr><tr><td>A</td><td>0.5</td><td>ζ</td><td>1.5</td></tr><tr><td>Ω</td><td>0.3</td><td>ζ<sub>p</sub></td><td>32.0</td></tr></tbody></table>

the spontaneous splay induced by motor activity.

V. ASTER FORMATION IN BULK ACTIVE POLAR MODEL

The gain more insight into the formation of multiple asters we study a simplified model for polar active fluids in bulk. We consider the following equations for motor concentration and polarization field,

$$ \partial_t m = D_m \Delta m + v_0 \nabla \cdot (m\mathbf{p}), \qquad (S13) $$

$$ \partial_t \mathbf{p} = -\zeta_p \nabla m + \Gamma_p (\kappa_p \Delta \mathbf{p} + \alpha_p \mathbf{p} - \beta_p |\mathbf{p}|^2). \qquad (S14) $$

Here, we use $\Delta \equiv \nabla^2$ is the Laplace-Beltrami operator In this simplified model, we consider a spontaneous long range order in the polarization instead of nematic order. We also neglect motor binding kinetics, and instead set the total number of motors to a constant average value, $m_0$, which corresponds to the steady state of motor binding kinetics: $m_0 = k_{on} \psi_0 / k_{off}$. Further, we suppress the density kinetics by going to the incompressible limit and setting the compressibility $\zeta_0 = 0$.

Equations (S13)-(S14) are integrated using a finite difference method with periodic boundary condition in a square box with side length $l$. For the time integration, we use a fourth order Runge-Kutta method with a timestep, $\delta t = 0.01$. We initialize the system with an aster in the polarization field $\mathbf{p} = (\cos\phi, \sin\phi)^T$ and a uniform motor concentration, $m_0$. All simulation parameters are given in Table IV.

First, we study the system size dependence by varying the box size $l$, and initial motor concentration $m_0 = 0.5$. In Fig. S3 the number of asters as a function of system size $l$ is shown. The single aster state from the initial condition is destabilized above a critical system size $l_c = 54$ and we find a monotonic increase in number of asters. Linearizing Eq.(S13)-(S14) around the polarized state $\mathbf{p} = \mathbf{e}_y + \delta\mathbf{p}$, $m = m_0 + \delta m$ and transforming into Fourier space with wavevector $\mathbf{k} = (k_x, k_y)^T$ gives

$$ \partial_t \delta \tilde{m} = -D_m |\mathbf{k}|^2 \delta \tilde{m} + i v_0 m_0 \mathbf{k} \cdot \delta \tilde{\mathbf{p}} + i v_0 i k_y \delta \tilde{m}, \qquad (S15) $$

$$ \partial_t \delta \tilde{\mathbf{p}} = -\zeta \mathbf{k} \delta \tilde{m} - \Gamma_p \kappa_p |\mathbf{k}|^2 \delta \tilde{\mathbf{p}}, \qquad (S16) $$
<table><caption>TABLE IV. Parameter values used in bulk simulations.</caption><thead><tr><th>Parameter</th><th>Value</th></tr></thead><tbody><tr><td>Δx Δy</td><td>1.0 1.0</td></tr><tr><td>D<sub>m</sub></td><td>0.5</td></tr><tr><td>v<sub>0</sub></td><td>0.35</td></tr><tr><td>ζ<sub>p</sub></td><td>10.0</td></tr><tr><td>κ<sub>p</sub></td><td>0.5</td></tr><tr><td>α<sub>p</sub></td><td>9.0</td></tr><tr><td>β<sub>p</sub></td><td>10.0</td></tr><tr><td>Γ<sub>p</sub></td><td>10.0</td></tr></tbody></table>

<figure><img src="image_7.png" /><figcaption>FIG. S3. Number of aster as a function of system size.</figcaption></figure>

where we neglected the Landau terms. This system has a maximal eigenvalue whose wavevector is given by

$$k_{\max} \approx \frac{\sqrt{m_0 v_0 \zeta}}{D_m + \Gamma_p \kappa_p}. \qquad (S17)$$

The corresponding wavelength $L_{max} = 2\pi/k_{max}$ determines the grid spacing of asters in our simulation. Hence, we expect to find more than one aster when the system is larger then $2L_{max} = 52.2$ which agrees reasonably well with the system size we find in our simulations (see Fig. S3).

Second, we investigate the number of motors in the system while keeping the system size constant at $l = 32$. We initialized the system with an aster in the polarization field, and varied the initial motor concentration $m_0$. Figure S4 shows the number of asters as a function of motor concentration. We find a critical motor concentration $m_0 = 1.5$ above which there are multiple asters emerging. Using the maximal wavevector determined Eq. (S17) by our linear stability analysis we can compute the motor concentration at which the corresponding wavelength $L_{max}$ is half the system size, such that it can support two asters. The resulting critical motor concentration,

$$m_0 = \frac{4\pi^2(D_m + \kappa_p \Gamma_p)^2}{\zeta(l/2)^2 v_0}, \qquad (S18)$$

has a value of $m_0 = 1.3$, which agrees reasonably well with the critical motor concentration we find from our simulations.

In Fig S5, we show the critical motor concentration at which one aster becomes unstable as a function of systems size. We show critical concentrations from both our simulations and our linear stability analysis, which agree well. The critical motor concentration becomes smaller with larger system size.
<figure><img src="image_8.png" /><figcaption>FIG. S4. Number of asters as a function of bound motors $m_0$ in a bulk simulation.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>FIG. S5. Critical motor concentration as a function of system size in a bulk simulation. Circles show the critical motor concentration from our simulations and the solid line show the result of our linear stability analysis.</figcaption></figure>

## VI. STEADY STATE ASTER SOLUTION

In order to find a steady state solution corresponding to an aster, we focus on the motor and polarization equations, which at steady state are given by,

$$0 = -\zeta_p \nabla m + \Gamma_p \kappa_p \Delta \mathbf{p}, \qquad (S19)$$

$$0 = \nabla \cdot [D_m \nabla m + v_0(m\mathbf{p})]. \qquad (S20)$$

Here, we neglected binding kinetics in the motor equation as well as the decaying terms for the polarization $\mathbf{p}$ and set compressibility $\zeta_0 \to 0$. A solution to Equation (S20) can be obtained by satisfying, $\nabla \ln(m) = -v_0/D_m \mathbf{p}$. We insert this into the Laplacian of Eq.(S19) and find a solution,

$$0 = \lambda_m m + \Delta \ln(m), \qquad (S21)$$

where $\lambda_m = \frac{\zeta_p v_0}{\Gamma_p \kappa_p D_m}$. This is known as the Liouville-Bratu-Gelfand equation [50-52] and has been solved for a number of physically relevant cases [53]. A single aster solution is given by

$$m = \frac{8}{\lambda_m (|\mathbf{r}|^2 + \xi_c^2)^2}, \qquad (S22)$$

$$\mathbf{p} = p_0 \frac{4\mathbf{r}}{|\mathbf{r}|^2 + \xi_c^2}, \qquad (S23)$$
with $\xi_c$ as the characteristic aster size and $p_0 = \frac{D_m}{v_0}$.

This set of solutions allows for $mp = \nabla\omega$ with

$$ \omega = -\frac{8p_0}{\lambda_m(\xi_c^2 + \mathbf{r}^2)^2} \qquad (S24) $$

We can then write the active term in the density kinetics of Eq. (5) as $-\zeta\nabla \cdot (mp) = -\zeta\nabla^2\omega$. The density kinetics are then expressible by the variational derivative of an effective free energy, $F_{\text{eff}} = -(\zeta/\Gamma_{\psi}) \int d^2x \omega\psi + F_{\psi}$.

Using the expression (S24) for $\omega$ in the free energy term, we can obtain the contribution of the $\zeta$-dependent active flux term to the droplet surface energy. By assuming that $F_{\psi}$ determines the circular droplet shape, such that $\psi(\mathbf{x}) = 1$ for $|\mathbf{x}| < R$ and zero otherwise, we can calculate this contribution as,

$$ \int d^2x \omega\psi = 2\pi \int_0^R dr \frac{8p_0}{\lambda_m} \frac{rdr}{(\xi_c^2 + \mathbf{r}^2)^2} = \frac{8\pi p_0}{\lambda_m \xi_c^2} \frac{R^2}{(\xi_c^2 + R^2)}, \qquad (S25) $$

which in the limit of small droplet size compared to aster size, gives a total contribution to the effective free energy, $-(\zeta/\Gamma_{\psi}) \int d^2x \omega\psi \simeq -\zeta \cdot \frac{8\pi p_0}{\Gamma_{\psi}\xi_c^4} R^2$. Including the equilibrium surface energy of a circular droplet, we get for the total effective free energy that determines droplet density, $(\gamma - \zeta \cdot \frac{8\pi p_0}{\Gamma_{\psi}\xi_c^4})R^2$, which leads to droplet instability when $\zeta$ is large enough.