Soft channel formation and symmetry breaking in exotic active emulsions

We use computer simulations to study the morphology and rheological properties of a bidimensional emulsion resulting from a mixture of a passive isotropic fluid and an active contractile polar gel, in the presence of a surfactant that favours the emulsification of the two phases. By varying the intensity of the contractile activity and of an externally imposed shear flow, we find three possible morphologies. For low shear rates, a simple lamellar state is obtained. For intermediate activity and shear rate, an asymmetric state emerges, which is characterized by shear and concentration banding at the polar/isotropic interface. A further increment in the active forcing leads to the self-assembly of a soft channel where an isotropic fluid flows between two layers of active material. We characterize the stability of this state by performing a dynamical test varying the intensity of the active forcing and shear rate. Finally, we address the rheological properties of the system by measuring the effective shear viscosity, finding that this increases as active forcing is increased—so that the fluid thickens with activity.


Scientific RepoRtS
| (2020) 10:15936 | https://doi.org/10.1038/s41598-020-72742-9 www.nature.com/scientificreports/ contractile swimmers (or pullers) behave as in-warding force dipoles and the flow structure is reversed. Analytical and numerical works confirmed that extensile suspensions can strengthen the imposed flow, thus leading to the reduction of the viscosity of the fluid, while contractile systems have the opposite effect and they increase the viscosity, a phenomenon that is commonly addressed as active thickening 28 . Most experiments concerning the rheology of active fluids focused on the behaviour of a single fluid system, where the active constituents are uniformly distributed. Nevertheless, many experiments have shown that it is now possible to encapsulate active material in a water-in-oil emulsion 24,25 to have droplets or shells of active material. Confining active liquid crystals is key to control active fluids and to localize the typical scales of energy injection. Moreover, it constitutes by itself a first step towards the development of micro-cargos for drug-delivery.
Simulations based on active gel theory have been proved to be able to qualitatively replicate the dynamics of active droplets-from the regular motion of topological defects 8,29 , to the turbulent regime 15,30 -and predict new motility modes 29,[31][32][33] . In this paper we consider a theoretical model of an active binary mixture [34][35][36] , to numerically study the rheological response of contractile mixture under an externally imposed shear.
The emulsion is prepared by mixing an isotropic passive fluid with a polar active one in presence of a surfactant which favours the emulsification of the two phases. We varied systematically the intensity of the active forcing and the shear rate, integrating the equations by means of a hybrid lattice Boltzmann approach 35,37,38 . We found that competition between contractile activity and the imposed shear leads to interesting phenomena: in particular weak contractile activity is able to order the lamellar pattern when shear rate is absent or highly suppressed, while for more intense active forcing the system undergoes symmetry-breaking in the direction normal to the walls leading to pronounced concentration banding-a phenomenon that we also observed in the passive counterpart at high shear rates. Interestingly, if both external and active forcing are strong enough, the symmetry is restored and a channel of passive isotropic fluid is formed in the middle of the system, while two soft active layers adhere to the wall. We tested the stability of such configuration by varying both activity and shear rate over time, finding that contractile forcing is able to stabilize the channel, even starting from an asymmetric state. Finally, we found that the effective viscosity of the system increases together with the intensity of the contractile activity. Such dynamic behavior is drastically different from that of the extensile mixture under shear, which is found to exhibit a multi-stable dynamics characterized by states at null and negative viscosity 36 , in agreement with the experimental results of Lopez and Guo.
The article is organized as follows. In the next section we will present the model used to describe our system (further details on the numerical methods and simulation parameters can be find in the "Appendix"). Later on, we will provide a morphological characterization of the various states obtained at varying activity and shear rate and we will thoroughly address the dynamics of the symmetry breaking and the soft channel formation. Finally we will present the results of the dynamical study performed varying activity and shear rate over time. The rheological measurements are outlined in the final discussion.

the model
In this Section we will present the physical variables relevant for the characterization of the nemato-hydrodynamic state of the system, its equilibrium properties in absence of activity and the evolution equations ruling the dynamics. We introduce the scalar field φ(r, t) as the concentration of active material and the polarization P(r, t) -a vector field accounting for the mean orientation of the active constituents. The velocity of the fluid will be denoted by v(r, t) and the density by ρ(r, t).
The equilibrium properties of the system are encoded in a modified Landau-Brazovskii 34-36,39,40 free-energy functional F [φ, P] that can be expressed as the first contribution defines the thermodynamic property of the binary mixture: The first term allows for the segregation of the two phases when a > 0 , as the free energy has two minima at φ = 0 and φ 0 , respectively corresponding to the passive and active phase, while φ cr = φ 0 /2 . The gradient terms determine the surface tension. Here we choose k φ < 0 so to favour the formation of interfaces throughout the system, while c must be positive to guarantee thermodynamic stability. This is equivalent to dispersing a suitable amount of surfactant with a relaxation time much smaller than either of the other two phases. For a symmetric composition of the mixture the theory defined by F bm [φ] admits as ground state a lamellar modulation with wave-number κ = |k φ |/2c when a < k 2 φ /4c . The second contribution to F is borrowed from the theory for liquid crystals and adapted for the treatment of a vector order parameter 41 Here the first two terms define the bulk properties of the liquid crystal that is confined in those regions where φ > φ cr , while the term proportional to k P describes the energy cost due to elastic deformation of the polar pattern. Finally, the last term on the right hand side of Eq. (1) defines the anchoring interaction between the concentration and the polarization field www.nature.com/scientificreports/ This contribution favours homeotropic anchoring of P at the interfaces. In particular, when β > 0 , the polarization points towards the decreasing values of φ . The ground state defined by the free energy functional F [φ, P] is characterized by a transition from the isotropic towards the lamellar phase for a < k 2 φ /4c + β 2 /k P 36 , with the polarization P confined in only one of the two phase and normally anchored to the interfaces. We note that the presence of a strong homeotropic interface anchoring is crucial to ensure the stability of the lamellar phase in active mixtures. Indeed, if β = 0 , lamellae may get disrupted by the stirring effect of the fluid flow, since the underpinning effect of the polarization field disappears 34 .
The dynamic equations governing the physics of the system are in the limit of incompressible fluid (which is a good approximation for active fluids). Equation (5) is the Navier-Stokes equation, where P is the ideal pressure and the stress tensor has been divided in a passive and an active contribution. The former is due to dissipative (viscous) and reactive phenomena resulting from the dynamics of the two order parameters. Its explicit expression is given by the sum of three terms whose explicit definition is given in the "Appendix". The active contribution has phenomenological origin and it is the only term not stemming from the free energy. It accounts for the local stress exerted by the active particles and is given by 3,42 where ζ is the activity parameter, positive for extensile systems and negative for contractile ones. By assuming the amount of active material to be conserved, φ evolves according to a convection-diffusion equation (Eq. (6)) in which M is the mobility and µ = δF /δφ is the chemical potential.
The evolution of the polarization field is governed by the Ericksen-Leslie equation for a vector order parameter (Eq. (7)), where Ŵ is the rotational viscosity, while D and ˜ respectively represent the symmetric and the antisymmetric part of the velocity gradient tensor ∂ α v β .
Moreover, we observe that a further self advecting term of the form w∇ · (φP) may be included on the lefthand side of the evolution equation for the concentration field-with w a reference speed-together with an analogous term in the Ericksen-Leslie equation for P of the form wP · ∇P . Such terms are commonly used to incorporate the description of individual particles able to swim in the fluid, hence, to move with respect to the fluid with a typical speed w. In the following, we will only consider systems with w = 0 . This assumption corresponds to consider shaker particles, namely organisms only capable to stress the surrounding fluid.
The simulations were performed on a bi-dimensional lattice, initializing the concentration uniformly as φ = φ cr . The polarization was randomly initialized both in orientation and intensity. We confined the system in a channel of width L with no-slip boundary conditions at the top and the bottom walls ( y = 0 and y = L ). These conditions were implemented in the LB algorithm by bounce-back boundary conditions and periodic boundary conditions in the y direction. The flow is driven by moving walls, with velocity v w for the top wall and −v w for the bottom wall. This gives a shear rate γ = 2v w /L . Neutral wetting boundary conditions were enforced by requiring on the wall sites that: where ∇ ⊥ denotes the partial derivative computed normally to the walls and directed towards the bulk of the system. The first condition ensures density conservation, and the second determines the wetting to be neutral. Concerning the boundary conditions for the polarization we impose strong tangential boundary conditions where P ⊥ and P denote, respectively, the normal and tangential components of the polarization field with respect to the walls. This is because, according to experiments, both bacteria and cytoskeletal extracts tend to orient along the wall directions when they are close to the boundaries.
Importantly, we observe that, although an externally imposed shear-flow can give rise to a first-order transition from isotropic to polar through a renormalization of the coefficients for the double-well potential of the polarization field 43 , such scenario does not occur for our choice of parameters. This is essentially because the shear rate is never strong enough to induce a change in the sign of the coefficient of the quadratic term www.nature.com/scientificreports/ proportional to α , which remains positive in those regions where φ < φ cr and the transition does not occur. Other numerical details concerning the numerical method and the mapping to physical values are given in the "Appendix".

Results
Asymmetry and soft channel formation in contractile emulsions. In order to analyze the behavior of our system, we performed a systematic scanning of both contractile activity and shear rate. In this Section we will present the morphological properties, summarized in Fig. 1, showing some examples of the different regimes occurring in our system. At equilibrium ( ζ = 0 ) the free energy in Eq. (1) prescribes the formation of a lamellar configuration for the concentration field. In absence of forcing (either external or internal), and starting from a mixed state, large domains of well aligned lamellae develop throughout the system. Complete ordering is however unlikely to be reached due to persistent dislocations that diffusion or backflow are not able to flush away. If a weak shear rate is externally imposed ( γ 0.8 × 10 −4 ), the lamellar texture aligns in the flow direction as shown in panel (a) of Fig. 1. In this case, dislocations in the lamellar pattern still persist in the long time dynamics. This changes when an active emulsion is considered. By keeping fixed the shear rate and switching on weak contractile activity ( |ζ | 4 × 10 −3 ), a completely ordered configuration is achieved with no defects in the lamellar arrangement, regardless of the initial conditions (see Fig. 1b). This is due to the mutual interaction between the external shear and the activity: while the former has the effect of aligning the pattern to the flow, the latter provides a nonvanishing energy contribution which keeps on stirring the fluid and eventually leads to the elimination of defects.
A further increment of the intensity of the active forcing engenders an important morphological transition. When activity is strong enough, the lamellar texture becomes unstable to deformations of the polar pattern. The underlying mechanism is the generic instability of active contractile materials to splay distortions, which lead to the formation of a self-assembled emulsion of passive droplets in an active background. This is energetically favoured, with respect to the lamellar configuration, as it allows for relieving elastic stresses. More interestingly, this phenomenology-which was already observed in unconfined systems 34 -is accompanied by symmetry breaking along the gradient direction. While an active matrix interspersed with passive droplets gathers and grows on one side of the channel, the rest of the system is left in a state poorly populated of active material. This www.nature.com/scientificreports/ eventually leads to the pearlification of the lamellar pattern 44 which results, in the long term dynamics, into a layer of isotropic fluid, with only some active features scattered with no definite pattern (see panel (c) of Fig. 1 and Suppl. Movie 1). As we will discuss in more depth in the following Section, for low shear rate and intense activity, the symmetry breaking is the result of the advective transfer of active material from one side of the channel to the other.
To investigate the role of the external forcing in such dynamical process, we considered the case of stronger shear rates (see panels (d)-(f) in Fig. 1). We found that, even a passive emulsion can undergo a similar morphological transition to an asymmetric state, as suggested by the contour plot of the concentration field in panel (d). Nevertheless, the matrix of polar liquid crystal is now more uniform if compared to its active analogue, and droplets of isotropic fluid only develop far from the wall, where the polarization roughly aligns at the prescribed Leslie angle 43 . More specifically, under shear flow, a passive polar gel is tilted with respect to the flow at an angle . We checked that, in the passive limit, the polarization aligns at the Leslie angle when the symmetry is broken and a liquid crystal slab is formed at one of the walls. This suggests that, for strong enough shear rates, the flow-aligning nemato-hydrodynamic effect becomes dominant and overtakes the relaxation towards the lamellar pattern favored by the free-energy. This gives rise to large polar domains with a well-aligned liquid crystal.
We then considered the case of an active contractile emulsion under a strong shear rate ( γ = 7 × 10 −4 ). Remarkably, under these conditions, activity prevents the development of any asymmetry (Fig. 1e). Two active layers form at the walls, with the proliferation of passive droplets, while a channel of isotropic fluid, populated by small active domains, develops in the centre of the system (Suppl. Movie 2). As activity is increased (see panel (f)) the channel shrinks and gets rid of polar droplets, while the active layers widen, due to the absorption of dispersed active material floating in the isotropic fluid.
To summarize, the asymmetry may develop due to both activity and external shear separately, whose mutual non-trivial interaction drives the growth of the liquid crystal domain sustained by the advection of the concentration field. Nevertheless, when both active and external forcing are strong enough, the competition of the two mechanisms leads to a symmetric configuration (see also the phase diagram in Fig. 5a). In the following we will refer to the shear as weak, whenever the lamellar phase is preserved ( γ 0.8 × 10 −4 ); analogously strong shear corresponds to a situation where the formation of the soft channel is favoured ( γ 7 × 10 −4 ).
In the next sections we will first thoroughly address the dynamical mechanisms driving the instability in the case of an active emulsion at strong activity and shear rate. Later on we will discuss the formation of soft channels, then the stability of the asymmetric state. the dynamics of symmetry breaking. In the previous Section we commented on the development of symmetry breaking in the concentration pattern occurring for weak shear rates and strong enough contractile activity. Panels (a)-(d) of Fig. 2 and Suppl. Movie 1 show the different stages of the process for the case with ζ = −0.002 and γ = 6.2 × 10 −4 . Starting from a totally mixed configuration, the system rapidly relaxes into a lamellar pattern. Due to the imposed shear, lamellae grow increasingly thinner until they disrupt into a pearlification texture (panel (a)). Concurrently, the excess of active material expelled from the bulk is progressively captured by the thick active matrix adhering on the walls. Initially both layers progressively widen (panel (b)), even if the absorption process may proceed at different speed, thus initiating the symmetry breaking.
When the lamellar texture disappears from the centre of the system, an emulsion of small active droplets remains suspended in the bulk of isotropic fluid. These are advected by the flow, whose structure is determined by both the imposed shear and the distribution of active material. To fully analyse the flow structure, we plotted the reduced velocity field ṽ in panel (f)-namely the difference between the actual velocity v and the imposed shear flow V-which exhibits a large convection roll. At an initial stage, this large vortex advects the droplets of active materials in closed trajectories from the upper half of the channel towards the lower half (see Movie 1). Eventually, one of the two layers ends up absorbing some of the active domains floating in the isotropic channel, thus breaking the top-bottom symmetry in the velocity profiles (see Fig. 2e). Once the symmetry is broken and the flow profiles are no longer symmetric, the vortex moves from the centre of the channel to the inversion region ( �v y � = 0 ), finally leading to the fully asymmetric configuration of Fig. 2d. This process is characterized by a net flux of the concentration field (see Fig. 2g). The flux here is computed as the amount of active material that flows across the mid-line of the simulation box in a lapse of 10 4 iterations ( � = dtdyφv z ). We notice that in the early dynamics-corresponding to the stage at which the lamellar pattern undergoes pearlification leading to the formation of the isotropic channel-the concentration flux is approximately null, while it attains a positive value during the shrinking of the bottom layer, showing that convection is responsible for the transfer of active material from the lower to the upper half of the channel. Finally, the flux decreases until it drops back to 0, and the system settles in the final asymmetric configuration of panel (d), characterized by a channel of isotropic fluid in the lower half and an active layer in the upper part of the system.
To quantitatively address the rheological features of our system, we measured the concentration profiles φ and the velocity profiles v averaged along the channel associated with the different stages of the symmetry breaking dynamics (see panels (e) and (h)). Two shear bandings are clearly visible as long as the system is symmetric and are associated with corresponding concentration banding (red and yellow lines). We observe that the effective shear rate is considerably lower in the regions where the inverse active emulsion is well developed, while it increases significantly in the centre of the channel. This behavior is related to the active nature of the matrices adhering to the walls. Indeed, the active liquid crystal produces in these regions a spontaneous flow which adds to the externally imposed shear flow, eventually thwarting it. As the bottom active layer melts, the velocity profile gets smoother and smoother so that the shear banding in proximity of the bottom wall progressively disappears and completely vanishes when the system becomes fully asymmetric. www.nature.com/scientificreports/ The late time configuration is characterized by strong concentration and shear banding, visible in the blue profiles in panels (h) and (e). In particular, we observe that the local shear rate is greater in the isotropic domain than in the active region, where the slope of the blue velocity profiles is lower-thus indicating that contractile activity induces an effective viscosity increase, i.e. active thickening.

Soft channel formation and active thickening.
In the previous section we discussed how activity triggers the breaking of the top/bottom symmetry, due to the development of convection rolls-sustained by the energy injected by the active constituents-that advect the concentration field, leading to a net flux that moves material from one side of the channel to the other. This transfer mechanism is effective for strong active forcing at low shear rates, since the strength of the activity-induced flows is comparable with that of the external forcing. But what happens if we are to raise the shear rate?
In this case the externally imposed shear flow strengthens the top-bottom symmetry, thus leading to a situation where the symmetry breaking would not occur and the system sets into a configuration characterized by a channel of isotropic fluid confined between two layers of active material. Moreover increasing the contractile activity stabilizes the channel (see for instance the phase diagram in Fig. 5). Indeed, the morphological transition from the asymmetric towards the soft channel state takes place at progressively lower values of the shear rate γ as the intensity of contractile activity is raised.
To show the behavior of our system, we now consider the case where both shear and active forcing are strong. The system, initialized in the mixed configuration, undergoes a similar early-time dynamics as described in the previous section: first a lamellar pattern is formed, then the texture is progressively destroyed due to pearlification instabilities. Two layers, characterized by a suspensions of passive droplets in an active background, develop at the walls and progressively thicken due to absorption of the active material, resulting in a symmetric configuration (see contour plots in the inset of Fig. 3b and Suppl. Movie 2 for the case at γ = 6.3 × 10 −4 and ζ = −0.006 ), where a soft channel of isotropic fluid forms in the centre of the system with a few suspended active features. This configuration is found to be stable even at long simulation times, since the strong shear flow prevents convection to set up an effective concentration flux. Interestingly, the shear rate plays an important role in controlling the amount of active material dispersed in the isotropic channel. Indeed, as the external forcing is increased the number of small active droplets progressively lowers, in a similar fashion to what happens when contractile activity is increased (see panels (e) and (f) of Fig. 1). More quantitatively, we plotted in Fig. 3a the concentration profiles for two values of v w , by keeping fixed the activity ( ζ = −0.006 ). We observe that, in both cases, the φ profiles exhibit strong fluctuations at the boundaries, due to the dynamical nucleation of the passive droplets via a mechanism resembling Ostwald ripening. In this region, active material is denser for the case at higher shear rate. As the interface between the active layer and the isotropic channel is approached, the profiles drastically  Fig. 3. Analogously to what happens in the case of asymmetric configurations discussed in the previous section, shear rate is more pronounced in the isotropic region, than in the active layers, where the effective viscosity is higher. However, this local thickening effect is also affecting the global rheological properties of the system, leading to the increase of the total effective viscosity with the intensity of the active forcing (see Fig. 5b).
We finally note that, if one considers a highly asymmetric mixture in which the active species are the minority components, the configuration here presented is inverted. Now a channel of active material flowing between two isotropic layers is obtained, with shear rate more pronounced close to the walls, while the shear-thickening feature is preserved.
This concludes the characterization of the dynamical regimes. A final question needs to be answered: does the system exhibit story dependence? Are there any hysteresis features? We will address this subject in the following section and we will show how the channel configuration is effectively stabilized by intense contractile activity even in presence of a partial asymmetry, given that this is not far developed.

Hysteresis features of soft channel formation.
In order to investigate the stability of the symmetric configurations we performed a dynamical study varying activity and shear rate over time. We first initialized the system in the mixed state and we started the simulation imposing a shear rate ( γ = 7.0 × 10 −4 ) such to obtain an asymmetric configuration. As the system develops an appreciable asymmetry (as shown in the top left panel of Fig. 4) we turned on the activity ( ζ = −0.007 ). We chose a value corresponding to a case characterized by the formation of an isotropic channel and we waited for the system to settle into a stable configuration. The final state is symmetric (central top panel), proving that a strong contractile activity is able to stabilize the channel, even starting from an asymmetric state. Importantly, this configuration is sustainable only if an external shear flow is imposed. Indeed, switching off the external forcing leads to the lamellar state (top right panel of Fig. 4) that one would obtain if the system were initialized in the mixed state. The starting asymmetric configuration can be now obtained by switching the activity off and the shear rate on at its original value. This dynamical test shows that contractile activity is fundamental to stabilize the isotropic channel, even if an applied external forcing is necessary for this configuration to be stable. Nevertheless, our system still exhibits some history dependence. Indeed, in the first stage of our test, we turned on the activity when the asymmetry was not fully developed. What would happen if we perform the same test starting from a more asymmetric state, as the one shown in the bottom central panel of Fig. 4? In this case the amount of active material in the upper part of the system is not enough to set up an effective counteracting flow and the system develops an asymmetric configuration, analogous to those observed at lower shear rate (bottom right panel). Nevertheless, by switching off the shear rate once more, it is possible to recover the lamellar state.

Discussion and conclusions
In this work we analysed the morphology and rheology of a mixture of an isotropic fluid and a polar active one under shear. We implicitly modelled the action of a surfactant which favours emulsification in the absence of shear. By varying the strength of the activity parameter (which measures the force generated internally by the www.nature.com/scientificreports/ system) and the imposed shear rate, we showed that the system can attain three possible morphologies (see phase diagram in Fig. 5a). For low shear rates ( γ 0.8 × 10 −4 ), the emulsion settles into a lamellar pattern. Defects in this regime are eliminated by increasing activity, giving rise, for sufficiently strong activity, to a state where lamellae align along the direction of the channel. For larger shear rates, we observe a morphological transition to a completely different state, characterised by phase separation where an active matrix pierced by a passive droplet emulsion coexists with an isotropic fluid with active droplets. This phase separation is driven by a convective flux which relocates active material from one side of the system to the other; this regime is additionally characterised by a strong concentration and shear banding as the two phases flow differently and have quite distinct average composition. There are two morphologies corresponding to the phase separated state. At low activity, the structure is asymmetric, and the isotropic phase and the active matrix migrate to different sides.  www.nature.com/scientificreports/ For large activity, the structure is remodelled to create an internal isotropic channel flanked by two soft active matrices wetting the walls. Our simulations suggest that the formation of a symmetric channel requires activity, to dissolve the active droplets in the middle of the isotropic phase.
To characterise the rheological properties of the system, we plotted in Fig. 5 the behaviour of the apparent viscosity η , normalised by the background fluid viscosity (i.e., the one entering the purely viscous part of the stress tensor). At constant activity, the viscosity lowers when shear rate is increased, so the system is shear thinning. Additionally, we observe that activity leads to an increase in viscosity, as in single-phase contractile fluids 11 , an effect which we refer to as active thickening. Active thickening is more visible for small shear rates, while it becomes negligible when the wall velocity is very strong, as the externally imposed flow dwarfs the active one.
Finally, we observe that our results are expected to remain valid also if an active liquid crystal with nematic, rather than polar, symmetry is considered-being the kind of allowed topological defects the most relevant difference between the two. This is because most of the observed behaviors are due to the interplay between the elastic properties of the liquid crystal and the structures developed by the concentration field, whose patterning is mostly ruled, for the range of ζ here considered, by the typical Brazovskii length-scale which controls the width of the lamellar texture and eventually the size of the droplets of passive/active material. Conversely, the presence of a self-advection term in our equations (of the form w∇ · (φP) mimicking, for instance, self-propelled bacteria) is going to alter the morphology and, likely, the rheological features of the mixture. Indeed, we expect our results to stay valid as long as w ≪ v w and w ≪ 2πζ /(κη) , corresponding to the case where the swimming speed w is small compared to the wall speed and to the active flows generated by the active stress tensor in the Navier-Stokes equation, while we expect our system to end up in an isotropic mixed state in the limit of strong w.
Our study shows that activity and an imposed flow can be harnessed to enhance the self-assembly potential of active systems, and specifically to create soft microfluidic channels by simply mixing suitable ingredients. It would be of interest to perform fully 3D simulations of this system to see how spatial dimensionality affects pattern formation, and whether a soft channel can be self-assembled in this case as well. Additionally, it would be desirable to study active self-assembly in active mixtures and emulsions experimentally. A direct replication of the ideas proposed here may be challenging, as it would require the creation of a two-dimensional emulsion. However, similar studies may be possible with actomyosin-water emulsions in the presence of a surfactant. Our results suggests that systems like these are likely to display a very rich and intriguing phenomenology.

numerical details and mapping to physical values
The stress tensor appearing at the right-hand side of the Navier-Stokes equation has been divided in two parts, respectively addressed as passive σ pass and active σ act . The first one is in turn the sum of dissipative and reactive contributions. The viscous term is given by where η is the shear viscosity. The second term in Eq. (8) is an interface contribution and is given by where f is the free energy density. The last term gives the stress due to deformations of the liquid crystal pattern. Its explicit form is where h = −δF /δP is the molecular field 41 and the shape factor ξ controls the anisotropy of the liquid crystal molecules and its response to a shear flow. In particular, it selects rod-like particles if positive or disk-like ones if negative. Moreover, the polarization field exhibits flow aligning or flow thumbling features under shear if |ξ | > 1 or |ξ | < 1 , respectively.
Equations (5)-(7), in the limit of incompressible flow, are numerically solved by means of a well validated hybrid lattice Boltzmann method. The Navier-Stokes equation was solved through a predictor-corrector LB scheme. The evolution equations for the order parameters φ and P were integrated through a predictor-corrector finite-difference algorithm implementing first-order upwind scheme. Fourth order accurate stencils for space derivatives were used. Simulations are performed on a two-dimensional square lattice (D2Q9) whose linear size ranges from L = 256 to L = 512 . The system is initialized in a mixed state, with φ uniformly distributed between 1.1 and 0.9. The concentration φ ranges from φ = 0 (passive phase) to φ ≃ 2 (active phase), that correspond to the two minima of the double-well potential. The starting polarization field P is randomly distributed between 0 (passive phase) and 1 (active phase). Unless otherwise stated, parameter values are a = 4 × 10 −3 , k φ = −6 × 10 −3 , c = 10 −2 , α = 10 −3 , k P = 10 −2 , Ŵ = 1 , ξ = 1.1 , φ 0 = 2 , β = 10 −2 , η = 1.67 . By following previous studies 45 , an approximate relation between simulation units and physical ones (such as those of a contractile active gel) can be obtained by using as length-scale, time-scale and force-scale respectively the values L = 1 µ m, τ = 10 ms and F = 1000 nN (see Table 1). Throughout our simulations the Reynolds number, for the case in which the droplets are observed, is evaluated in terms of the average droplet radius, of the viscosity and of the velocity of the fluid. It remains below 0.25, a value in which inertial effects are indeed negligible. (P α h β − P β h α ) − ξ 2 (P α h β + P β h α ) − κ∂ α P γ ∂ β P γ ,