Active sorting to boundaries in active nematic–passive isotropic fluid mixtures

Abstract

We use a two-fluid model to study a confined mixture of an active nematic fluid and a passive isotropic fluid. We find that an extensile active fluid preferentially accumulates at a boundary if the anchoring is planar, whereas its boundary concentration decreases for homeotropic anchoring. These tendencies are reversed if the active fluid is contractile. We argue that the sorting results from gradients in the nematic order, and show that the behaviour can be driven by either imposed boundary anchoring or spontaneous anchoring induced by active flows. Our results can be tested by experiments on microtubule-kinesin motor networks, and may be relevant to sorting to the boundary in cell colonies or cancer spheroids.

---

1 Introduction

A central challenge in biophysics is to understand how numerous individual components self-organize to form more ordered structures. For example, experiments with starfish embryos in which endodermal and ectodermal cells were mixed showed that the cells spontaneously sorted, with endodermal cells migrating towards the interior and ectodermal cells moving to the outer boundary of the embryo.¹ Similarly, in Huang et al.,² a co-culture of homogeneously mixed metastatic breast cancer cells and non-tumorigenic epithelial cells initially sorted and later exhibited an inversion of tissue architecture. In confluent cell layers, different cell types can self-arrange into distinct regions,^3,4 and bacterial species within a biofilm have been observed to spontaneously sort based on differences in their mechanical properties.⁵

Biological phase separation has often been examined using the framework of equilibrium thermodynamics, where sorting is assumed to stem from the minimization of a free energy. Among the most prominent of the equilibrium theories is the differential adhesion hypothesis (DAH),^6,7 which attributes the ordering to variations in the binding strength between like and unlike cells. However, since living matter operates out of equilibrium, the DAH may not adequately account for factors such as chemical signalling, cell division, and mechanical forces. Indeed, recent experiments have demonstrated that the DAH alone is insufficient to fully explain cell sorting in living systems.^8–10

A class of systems exhibiting out-of-equilibrium phase separation is active matter,^11,12 which describes self-motile particles. A paradigmatic example of active phase separation is motility-induced phase separation (MIPS), where a system of self-propelled particles separates into dense and dilute phases.¹³ Other routes to active phase separation include aligning torques, chemical bonds, motility differences, and hydrodynamic interactions.^14–20 Active phase separation has been documented in various continuum models,^21–24 as well as in phase field and vertex model simulations of heterogeneous cell layers.^4,25,26 The sorting can also be influenced by confining the active material; in many active mixtures, the different species preferentially aggregate at or move away from a boundary.^27–29 Swimming microorganisms, in particular, congregate at walls by hydrodynamic interactions.^30,31

A widely used formalism for studying dense cellular aggregates is that of active nematics.^32–34 These comprise elongated particles that align in a common direction but lack positional order. The particles exert local dipolar forces on their surroundings by pushing or pulling along their long axis. A consequence of this is that active nematics in bulk are unstable to a chaotic flow state called active turbulence.³⁵ By contrast, in confinement, they can exhibit a rich variety of steady flow states, including unidirectional and bidirectional laminar flow, spontaneous rotational flows, and a lattice of vortices which can be controlled by adjusting the dimensions and the boundary conditions of the confining box.^36–39 Moreover, active nematics at a boundary or interface give rise to interesting new physics, notably active anchoring,⁴⁰ active waves,^41–43 finger formation,⁴⁴ and interfacial self-folding.⁴⁵

We have recently demonstrated that a binary mixture of an active nematic and a passive fluid will sort even in the absence of thermodynamic driving forces. However, many of the experiments demonstrating sorting are carried out using confined cell collectives, and the effect of boundaries on the active sorting remains unclear. In this paper we use a two-fluid model to investigate this problem numerically.^23,46,47 We show that, in the absence of thermodynamically-controlled wetting, whether the active or the passive component preferentially sorts to the boundary depends on the director anchoring, and we give analytical arguments to explain this result.

The paper is structured as follows: in Section 2, we review the model introduced in Bhattacharyya et al.^46,48 to study a mixture of an active nematic and a passive isotropic fluid and give details of the numerical approach. In Section 3, we confine the mixture to a square box, imposing either homeotropic or parallel anchoring on the nematic director. We find that an extensile (contractile) active nematic component preferentially sorts to the wall if the anchoring is planar (homeotropic). We then show that the behaviour is similar if the wall anchoring itself results from activity (Section 4). In Section 5 we turn to confinement in a circle. This allows us to stabilise a circulating flow pattern and hence to compare the effects of curvature forces to those resulting from gradients in the magnitude of the nematic ordering. The final Section 6, summarises the results and suggests a possible mechanism for the inversion of tissue architecture observed in Huang et al.²

2 Model

We consider a mixture of two fluid components i = 1, 2 where i = 1 is an active nematic and i = 2 is a passive, isotropic fluid. Each component satisfies mass and momentum conservation:

∂_tρⁱ + ∇·(ρⁱuⁱ) = 0, (1)

∂_t(ρⁱuⁱ) + ∇·(ρⁱuⁱuⁱ) = f^visc,i + f^thermo,i + f^body,i + γϕ(1 − ϕ)(u³⁻ⁱ − uⁱ) (2)

where ρⁱ and uⁱ are the density and velocity of component i and ϕ is the concentration of component 1. The right hand side of eqn (2) describes the forces acting per unit volume of the fluid, and accounts for viscous (f^visc,i), thermodynamic (f^thermo,i), and local body forces (f^body,i). In addition, the term γϕ(1 − ϕ)(u³⁻ⁱ − uⁱ) models relative drag between the two fluid components, the strength of which is controlled by the coefficient γ.

It is convenient to reformulate the equations of motion in terms of a centre-of-mass (COM) fluid, and a relative fluid. The density field of the combined COM fluid is ρ^c = ρ¹ + ρ², and the COM flow field is u^c = ϕu¹ + (1 − ϕ)u². We also define a relative flow field δu = u¹ − u². We will work in the limit where the relative drag is the fastest mode of relaxation in the system i.e. γ is large. Thus, the relative flow is much smaller than the combined velocity of the fluid i.e. |δu| ≪ |u^c|. Adding eqn (1) and (2) for each fluid, and dropping terms of order (δu)², gives the equations of motion for the COM fluid:^46,48–50

∂_tρ^c + ∇·(ρ^cu^c) = 0, (3)

We now describe the terms on the right-hand side of eqn (4). The viscous term can be written as the divergence of a viscous stress tensor,

where νⁱ is the kinematic viscosity of fluid i. To simplify things further, we assume that ν¹ = ν² = ν. Since |δu| ≪ |u^c|, the viscous term for the combined fluid is then^49,50

The second term on the right-hand side of eqn (4) is the thermodynamic force derived from a free energy density [scr F, script letter F], which we will describe below. The thermodynamic force acting on each fluid component can be written as

f^thermo,i = −ρⁱ∇μⁱ, (7)

where μⁱ = δ[scr F, script letter F]/δρⁱ is the chemical potential of species i. The total thermodynamic force acting on the COM fluid can be written as the divergence of a stress tensor, σ^thermo,:^46,49

The final term on the right-hand side of eqn (4) is a body force that models the active and passive forces arising from the nematic nature of the first fluid component. The nematic order is described by a rank-2 nematic tensor field⁵¹

where S is the magnitude and n is the alignment axis of the nematic ordering. The dynamics of Q is governed by the Beris–Edwards equation^32,52

(∂_t + u¹·∇)Q + W = ΓH, (10)

whereis the co-rotation term describing the response of the orientation field to the traceless⁵³ strain rate tensor Ẽ¹ and vorticity tensor Ω¹ in the fluid. The flow-tumbling parameter λ describes the relative importance of strain and vorticity. Finally,is a molecular field describing relaxation to the minimum of the free energy at a rate depending on Γ, the rotational diffusivity.

Distortions in the nematic field exert a restoring stress on the underlying fluid, known as the elastic backflow term. This is given by³²

Importantly, the active nematic species exerts dipolar forces along the local elongation axis. This can be modelled by an active stress^32,35

Π^act = −ζϕQ, (14)

where the activity parameter ζ describes the strength of the dipolar forces and the sign of ζ determines their direction. Positive (negative) ζ corresponds to extensile (contractile) activity, where the elementary dipoles exert a net force outwards (inwards) along their long axis. Combining eqn (13) and (14), the total body force appearing in eqn (4) is

f^body,i = δ_i,1∇·(Π^backflow + Π^act). (15)

It now remains to define the free energy [scr F, script letter F] which is used to calculate the chemical potential and the molecular field. We choose⁴⁸

where the first term defines an isothermal equation of state for the combined fluid. Thus, the combined fluid is essentially incompressible, and density variations are resolved much faster than the timescale of system dynamics. However, the individual component fluids are compressible, allowing sorting.

The second term in eqn (16), in curly brackets, is the Landau free energy. The material parameters a and b describe the equilibrium state of the system in the absence of activity. If a > 0 and b > 0, equilibrium corresponds to the active nematic and passive isotropic fluids being homogeneously mixed.

The final term in eqn (16), in square brackets, is the Landau–de Gennes free energy describing the thermodynamics of a passive nematic. C and K are bulk material properties of the nematic, and S₀ is the magnitude of the nematic order in equilibrium. We take S₀ = 0 corresponding to a paranematic, i.e. no nematic order in the absence of activity or boundaries.

We solve eqn (3) and (4) for the combined fluid using an incompressible lattice-Boltzmann method.⁵⁴ Eqn (1) and (2) for the compressible active component are solved using a compressible lattice Boltzmann method.^46,49 Eqn (10) for the Q-field is integrated using a finite difference approach. More details of the simulation techniques can be found in Bhattacharyya et al.^46,48 where we used the algorithm to study the ordering of bulk mixtures of an active and a passive fluid, and mixtures of active fluids.

In this paper, we extend our previous work by confining an active–passive mixture. We consider two different geometries. In the first of these, the mixture is confined inside a square box of side L with solid walls. We impose free-slip boundary conditions for each fluid component at the box walls. This is simulated numerically by using a lattice Boltzmann bounce-back scheme. The orientation field is constrained to be perfectly anchored at the walls in either a planar or homeotropic configuration.

In the second case, the mixture is confined inside a fixed region by a third highly viscous fluid, with no-slip and no-penetration conditions at the boundary. The nematic field at the boundary is free, and its orientation is controlled by active flows. Numerically, we use a fixed phase field Ψ to distinguish between the inside (Ψ = 1) and the outside (Ψ = 0) of the confinement. In order to have a smooth boundary, Ψ is chosen to minimize the free energy F_Ψ = A_ΨΨ²(1 − Ψ)² + K_Ψ(∇Ψ)² with A_Ψ = 0.4 and K_Ψ = 0.1. The external confining fluid is chosen to have the same density as the fluid mixture inside, but it is highly viscous (ν^out = 10/6) and has strong frictional damping f^fric,i_Ψ=0 = −10uⁱ, ensuring no flow outside the confinement. At the boundary, Ψ varies smoothly from 1 to 0 over 4 lattice sites. The activity is localized inside the confined region by setting it to 0 if Ψ < 0.90. To deal with the finite width of the interface, we linearly interpolate the kinematic viscosity and friction coefficient for Ψ ∈ (0, 1). The velocity field is continuous at the interface, and zero outside the boundary.

We use the parameters ρ^C = 40, ν₁ = ν₂ = 1/6, γ = 0.1, a = 0.0001, b = 0, Γ = 0.33, C = 0.001, K = 0.02, S₀ = 0, λ = 0.7, and ζ ∈ [0.001, 0.012]. The system is initialized in a uniformly mixed configuration, with equal concentrations of active and passive material (ϕ = 1/2) everywhere. We initialize the nematic tensor in a random configuration with noise, and then simulate for 500 burn-in steps with no activity. Simulations are then run for at least 6 × 10⁵ timesteps, with snapshots taken every 1000 timesteps, until a steady state is reached. We quantify the typical concentration and fluctuations of the active fluid near the boundary, ϕ_edge, by measuring the mean and standard deviation of ϕ in a small region near the boundary, across different time snapshots. For a circular (square) confinement, these regions are given by annuli (square annuli) of width 2 units.

3 Active sorting to boundaries driven by imposed anchoring

We first study an equal mixture of an active nematic and a passive isotropic fluid confined inside a square box (L = 120). The velocity field satisfies a free-slip boundary condition, while the nematic director field at the walls is anchored in a planar (parallel to wall) or homeotropic (perpendicular to wall) configuration.

In the bulk of the box we observe active turbulence and microphase separation as detailed in Bhattacharyya et al.⁴⁶ Note, however that the snapshots in Fig. 1 and Movies S1 and S2, which are for a mixture with extensile activity ζ = 0.01, suggest that the active nematic component preferentially accumulates at the boundary if the anchoring is planar (panel a) but has a lower concentration at the boundary for homeotropic anchoring (panel b).

Fig. 1 Active nematic–passive mixture with imposed strong anchoring: (a) extensile activity, planar anchoring: for extensile activity (ζ = 0.01), with planar anchoring, the concentration of active nematic at the boundary wall is enhanced. The colourbar shows the concentration of active nematic fluid ϕ. (b) Extensile activity, homeotropic anchoring: if the anchoring is changed to homeotropic, the concentration of active nematic at the boundary is depleted. (c) Concentration of active nematic near the wall, averaged over time, for varying activity. Error bars show the standard deviation.

These observations are quantified in Fig. 1(c) where the concentration of the active fluid component near the wall is plotted as a function of activity for the two anchoring configurations. This figure also allows a comparison of sorting for a contractile active component, where the opposite behaviour is seen, with the wall concentration being suppressed by planar, and enhanced by homeotropic, anchoring. This comparison indicates that the concentration ordering at the boundary is driven by active forces.

Indeed, the boundary sorting can be explained by noting the form of the active stress, eqn (14), which implies that gradients of nematic order or material concentration generate active forces. Consider a nematic field, with orientation n, at a boundary, with unit vectors m and l denoting the normal and tangential directions to the boundary, as shown in in Fig. 2(a). For a gradient of nematic order or concentration ∇(ϕS) that points along m, the component of the active force normal to the boundary is F^norm = ζ|∇(ϕS)|(2(m·n)² − 1)m.⁴⁰ Its direction depends on the anchoring angle of the nematogens at the boundary, for example for an extensile nematic with planar anchoring F^norm points in the same direction as the gradient, ∇(ϕS). The gradients of S and ϕ play a similar role as an interface between two immiscible fluids, and the resultant forces look similar to those described by Coelho et al.⁵⁵

Fig. 2 Forces at a boundary: schematic representation of (a) a nematic at a boundary showing the vectors defined in the text. Background colours show the magnitude of nematic order, with yellow denoting high nematic order and blue denoting low nematic order. The gradient of nematic order points from blue to yellow. (b) Flows due to an extensile (ζ > 0) active nematic. (c) Flows due to a contractile (ζ < 0) active nematic. (d) and (e) Nematic order gradient and active forces at an anchored boundary. Extensile nematics are ordered at the boundary, but disordered in the bulk (shown by curved black lines), resulting in a gradient in nematic order, S. The resultant active force f_act points (d) towards the boundary for planar anchoring, and (e) away from the boundary for homeotropic anchoring.

Thus the accumulation of the active fluid at the anchored boundary can be attributed to the nematic order gradient between the nematogens influenced by the boundary anchoring and the more disordered bulk. For an extensile nematic with planar anchoring, the order gradient sets up active flows towards the boundary (Fig. 2(b) and (d)). This increases the concentration of extensile nematic at the boundary, setting up a concentration gradient, ∇ϕ, which reinforces the order gradient. The balance between the active force and these restoring forces controls the magnitude of the sorting. The sorting is opposed by the chaotic flows of active turbulence and by the Landau free energy, which favours a mixed state. Fig. 2(e) contrasts the case of homeotropic anchoring, where the gradient-induced flows are away from the wall. The change in behaviour can easily be explained by noting the change in the orientation n. The flow direction would be reversed for a contractile nematic (Fig. 2(c)).

For extensile systems, the bulk is turbulent, and the velocity field in the bulk is higher than that at the edges (for both planar and homeotropic anchoring). However, the velocity field alignment near the edges is uniform, while it is chaotic in the bulk (Fig. A1 in the SI). Similar results are observed for a box with imposed anchoring and no-slip boundary conditions (see Fig. A2 in the SI).

4 Active sorting to boundaries driven by spontaneous active anchoring

We next investigate whether sorting to a boundary can be observed even in the absence of imposed anchoring. We consider a square box containing an active nematic and a passive isotropic fluid component in equal amounts, and delineated by a third highly viscous confining fluid. Boundary conditions on the flow are no-slip, and the nematic field at the box boundary is free, so that its orientation is controlled by active flows.

The results are summarised in Fig. 3 for extensile activity and a box length L = 160 which is sufficiently large that the bulk shows active turbulence. We quantify the degree of sorting by measuring the peak of the average concentration profile of the active fluid near the edge of the box. At low activities, ζ = 0.004, the snapshot in Fig. 3(a) (see also Movie S3) shows a slight accumulation of the active component at the boundary. For higher activities, ζ = 0.01 (Fig. 3(b) and Movie S4), there is a clear sorting of active fluid to the boundary. Fig. 3(c) demonstrates the increase of the boundary accumulation with the activity. While this is significant, it is weaker than that caused by strong thermodynamic anchoring.

Fig. 3 Active nematic–passive mixture with active anchoring: (a) extensile activity, active anchoring: at a low activity, ζ = 0.004, the active material has a slightly enhanced concentration at the boundary. (b) Extensile activity, active anchoring: at higher activity, ζ = 0.01, bulk flows are chaotic, but the active anchoring, and hence the sorting is stronger. The colourbar shows the concentration of active nematic fluid ϕ. (c) Average concentration of active nematic near the wall for varying activity. Error bars show the standard deviation.

The boundary sorting that occurs, even in the absence of imposed boundary conditions, is a result of active anchoring, the alignment of nematogens at a boundary caused by active flows.^40,56 Returning to the argument in Section 3 and Fig. 2, the component of the active force tangential to the boundary is F^tang = 2ζ|∇(Sϕ)|(m·n)(l·n)l.⁴⁰ If the director at the boundary is unconstrained, this force creates flows which tend to rotate the orientation field into the planar configuration for extensile fluids. The active anchoring at the surface then acts in the same way as thermodynamic anchoring to establish a gradient of nematic order normal to the surface which leads to preferential accumulation of the active component at the surface. (If λ > 0, there are no active flows for contractile activity as we are considering a system with no nematic order in the passive state and the shear-induced nematic ordering required for active turbulence in this limit requires λζ > 0.⁵⁷)

Similar results are observed for a channel with spontaneous active anchoring and no-slip boundary conditions (see Fig. A3 in the SI). Figures showing the spatial variation of the active concentration ϕ is presented in Fig. A4 of the SI.

5 Active sorting to boundaries driven by correlated director gradients

Returning to the expression (14) for the active stress, we have so far considered the active forces due to gradients in concentration or nematic order,

However, gradients in the direction of the nematic order parameter may also lead to flows that can push the active fluid component radially,

F^curv = 2ζSϕ∇·(nn). (18)

To investigate the relative contributions to sorting from gradients in the magnitude or direction of the order parameter, we consider a circular confining region. As predicted theoretically in Woodhouse et al.⁵⁸ and demonstrated experimentally in Opathalage et al.,³⁸ for a small confining radius active turbulence is replaced by azimuthal flows with the gradients associated with the curvature of the director field driving spontaneous rotation. This is shown in Fig. 4(a) (see also Movie S5) for R = 30 and an activity ζ = 0.005. Fig. 4(b) and (c) (see also Movies S6 and S7) show similar plots for R = 60 and R = 90, where the bulk shows chaotic active flows.

Fig. 4 Active nematic–passive mixture in circular confinement: extensile activity, active anchoring: concentration profile of the active nematic fluid, ϕ (colourbar) and director (black lines) for (a) R = 30, (b) R = 60, (c) R = 90 for ζ = 0.005. Comparison of the radial forces due to gradients in the nematic order parameter and gradients in the director field as a function of the radial co-ordinate, r, for (d) R = 30, (e) R = 60, (f) R = 90. (g) Average concentration of active nematic near the wall for varying radii of confinement.

The variation with R of the concentration enhancement near the edge of the circle is plotted in Fig. 4(g). This increases to R ∼ 60 and then decreases slightly as the director ordering is disrupted by the turbulent bulk flows. The no-slip boundary conditions imply that the maximum shear-induced nematic order is achieved slightly away from the edge of the confinement. Therefore we calculate the concentration enhancement by averaging over an annulus centred on the peak of the radially-averaged active concentration profile.

To assess the relative contributions of the nematic ordering and director curvature gradients to the radial force we assume azimuthal symmetry, S = S(r), ϕ = ϕ(r), and n = (cos Θ, sin Θ) where Θ = θ + Θ₀, with Θ₀ a constant tilting angle. It then follows from eqn (17) and (18) that⁴⁶

f_order·r̂ = −ζ cos(2Θ₀)∂_r(Sϕ), (19)

These contributions, which can be calculated numerically, are plotted in Fig. 4(d)–(f) for R = 30, 60, 90. We find that f_order and f_curv contribute approximately equally to setting up the concentration profile. The error bars on these curves indicate that the forces fluctuate much more strongly for larger R because of the active turbulence in the bulk. Note that f_order·r̂ becomes negative near the boundary – this is a consequence of the active anchoring being suppressed by the no-slip boundary conditions.

6 Discussion

To summarise, we have shown that a mixture of an active nematic fluid and a passive isotropic fluid spontaneously sorts at a boundary. The sorting mechanism depends on gradients in the nematic order which are set up by boundary anchoring which can be externally imposed or result from active boundary flows. An extensile (contractile) active nematic preferentially accumulates at a boundary where the anchoring is planar (homeotropic) but migrates away from a boundary with homeotropic (planar) anchoring. In the absence of any imposed anchoring, spontaneous anchoring leads to an active species accumulating at the boundaries of its confinement.

Balancing thermodynamic driving with active forces resulting from confinement provides a way to control the positioning of an active fluid component. An example is shown in Fig. 5 (see also Movie S8). A fluid component is initiated at the centre of a confining region, stabilised by thermodynamic phase separation (by choosing a < 0, b > 0 in eqn (16)). If the inner fluid becomes active it elongates into an elliptical shape, driven by the active instability,⁴⁰ and then spreads into a rotating ring at the boundary. Indeed, cells sorting to the boundary of a confined region or droplet is ubiquitous in biology, and has been observed in embryonic amphibian cells,⁵⁹ mouse embryoid cells,⁶⁰ regenerating hydra,⁶¹ reconstituted starfish embryo,¹ and cancerous tumours.^2,62 Our results highlight a possible role of active forces in driving the cellular organisation, which may be particularly relevant when cells organize opposite to the DAH prediction. For instance, architecture inversion in Huang et al.² would be consistent with active nematic stresses associated with cell division, while cell sorting outcomes in Moore et al.⁶⁰ may be driven by active nematic stresses generated by E-cadherin.⁴

Fig. 5 Active sorting: extensile activity, active anchoring: (a) stable initial configuration in a thermodynamically phase-separated system (a < 0, b > 0). (b) When the inner species becomes active, it elongates into an elliptical shape. (c) On contact with the boundaries, the active material spreads along the confinement. (d) The active material accumulates into a ring at the edge of the container, which spontaneously rotates. The colourbar shows the concentration of the active species, ϕ.

In future work, it would be interesting to study the effects of confinements with locally changing curvatures, and how this affects the boundary distribution of active material. It might be possible to design confinement shapes such that an active nematic material would spontaneously accumulate in a specific region. Another interesting avenue for exploration would be the evolution of an three-dimensional active droplet with a deformable shape which contains an active–passive mixture. The feedback loop between active flows deforming the shape of the droplet and boundary sorting changing the activity distribution may provide a route to the spontaneous organization of complex patterns.

Conflicts of interest

There are no conflicts to declare.

Data availability

Supplementary information: Videos of the active sorting shown in Fig. 1, 3, 4 and 5. See DOI: https://doi.org/10.1039/d5sm00636h.

The code for the simulations can be found at https://github.com/saraswatb/BoundarySorting.

Acknowledgements

We thank D. Cropper, I. Hadjifrangiskou, A. Mietke, J. Rozman, S. Thampi and K. Thijssen for useful discussions and feedback. This research was supported in part by grant no. NSF PHY-2309135 to the Kavli Institute for Theoretical Physics (KITP). SB acknowledges funding from the Rudolf Peierls Centre for Theoretical Physics, and the Crewe Graduate Award. JMY acknowledges support from the ERC Advanced Grant ActBio (funded as UKRI Frontier Research Grant EP/Y033981/1).

References

1. S. Suzuki, I. Omori, R. Kuraishi and H. Kaneko, Dev., Growth Differ., 2021, 63, 343–353.
2. Y. L. Huang, C. Shiau, C. Wu, J. E. Segall and M. Wu, Biophys. Rev. Lett., 2020, 15, 131–141.
3. S. F. G. Krens and C.-P. Heisenberg, Curr. Top. Dev. Biol., 2011, 95, 189–213.
4. L. Balasubramaniam, A. Doostmohammadi, T. B. Saw, G. H. N. S. Narayana, R. Mueller, T. Dang, M. Thomas, S. Gupta, S. Sonam, A. S. Yap, Y. Toyama, R.-M. Mège, J. M. Yeomans and B. Ladoux, Nat. Mater., 2021, 20, 1156–1166.
5. E. R. Oldewurtel, N. Kouzel, L. Dewenter, K. Henseler and B. Maier, eLife, 2015, 4, e10811.
6. M. S. Steinberg, J. Theor. Biol., 1975, 55, 431–443.
7. R. A. Foty and M. S. Steinberg, Dev. Biol., 2005, 278, 255–263.
8. J. D. Tse, R. Moore, Y. Meng, W. Tao, E. R. Smith and X.-X. Xu, BMC Dev. Biol., 2021, 21, 2.
9. P. Heine, J. Lippoldt, G. A. Reddy, P. Katira and J. A. Käs, New J. Phys., 2021, 23, 043034.
10. S. Pawlizak, A. W. Fritsch, S. Grosser, D. Ahrens, T. Thalheim, S. Riedel, T. R. Kießling, L. Ostwald, M. Zink, M. L. Manning and J. A. Käs, New J. Phys., 2015, 17, 083049.
11. M. C. Marchetti, J.-F. Joanny, S. Ramaswamy, T. B. Liverpool, J. Prost, M. Rao and R. A. Simha, Rev. Mod. Phys., 2013, 85, 1143.
12. M. E. Cates and C. Nardini, Rep. Prog. Phys., 2025, 88, 056601.
13. M. E. Cates and J. Tailleur, Annu. Rev. Condens. Matter Phys., 2015, 6, 219–244.
14. J. Zhang, R. Alert, J. Yan, N. S. Wingreen and S. Granick, Nat. Phys., 2021, 17, 961.
15. S. Gokhale, J. A. Li, A. Solon, J. Gore and N. Fakhri, Phys. Rev. E, 2022, 105, 054605.
16. S. Thutupalli, D. Geyer, R. Singh, R. Adhikari and H. A. Stone, Proc. Natl. Acad. Sci. U. S. A., 2018, 115, 5403–5408.
17. A. Furukawa, D. Marenduzzo and M. E. Cates, Phys. Rev. E:Stat., Nonlinear, Soft Matter Phys., 2014, 90, 022303.
18. N. J. Lauersdorf, E. Nazockdast and D. Klotsa, Soft Matter, 2025, 21, 6132–6143.
19. C. P. Beatrici and L. G. Brunnet, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2011, 84, 031927.
20. J. M. Belmonte, G. L. Thomas, L. G. Brunnet, R. M. C. de Almeida and H. Chaté, Phys. Rev. Lett., 2008, 100, 248702.
21. J. Stenhammar, A. Tiribocchi, R. J. Allen, D. Marenduzzo and M. E. Cates, Phys. Rev. Lett., 2013, 111, 145702.
22. R. Wittkowski, A. Tiribocchi, J. Stenhammar, R. J. Allen, D. Marenduzzo and M. E. Cates, Nat. Commun., 2014, 5, 4351.
23. C. A. Weber, C. H. Rycroft and L. Mahadevan, Phys. Rev. Lett., 2018, 120, 248003.
24. A. Tiribocchi, R. Wittkowski, D. Marenduzzo and M. E. Cates, Phys. Rev. Lett., 2015, 115, 188302.
25. J. N. Graham, G. M. Zhang and J. M. Yeomans, Soft Matter, 2024, 20, 2955.
26. J. Rozman and J. M. Yeomans, Phys. Rev. Lett., 2024, 133, 248401.
27. X. B. Yang, M. L. Manning and M. C. Marchetti, Soft Matter, 2014, 10, 6477–6484.
28. L. Neville, J. Eggers and T. B. Liverpool, Soft Matter, 2024, 20, 8395.
29. C. Pérez-González, R. Alert, C. Blanch-Mercader, M. Gómez-González, T. Kolodziej, E. Bazellieres, J. Casademunt and X. Trepat, Nat. Phys., 2019, 15, 79–88.
30. A. P. Berke, L. Turner, H. C. Berg and E. Lauga, Phys. Rev. Lett., 2008, 101, 038102.
31. S. E. Spagnolie and E. Lauga, J. Fluid Mech., 2012, 700, 105–147.
32. A. Doostmohammadi, J. Ignés-Mullol, J. M. Yeomans and F. Sagués, Nat. Commun., 2018, 9, 3246.
33. T. B. Saw, A. Doostmohammadi, V. Nier, L. Kocgozlu, S. Thampi, Y. Toyama, P. Marcq, C. T. Lim, J. M. Yeomans and B. Ladoux, Nature, 2017, 544, 212–216.
34. G. Duclos, C. Blanch-Mercader, V. Yashunsky, G. Salbreux, J.-F. Joanny, J. Prost and P. Silberzan, Nat. Phys., 2018, 14, 728–732.
35. R. A. Simha and S. Ramaswamy, Phys. Rev. Lett., 2002, 89, 058101.
36. R. Voituriez, J.-F. Joanny and J. Prost, Europhys. Lett., 2005, 70, 404.
37. T. N. Shendruk, A. Doostmohammadi, K. Thijssen and J. M. Yeomans, Soft Matter, 2017, 13, 3853–3862.
38. A. Opathalage, M. M. Norton, M. P. N. Juniper, B. Langeslay, S. A. Aghvami, S. Fraden and Z. Dogic, Proc. Natl. Acad. Sci. U. S. A., 2019, 116, 4788–4797.
39. F. Caballero, Z. You and M. C. Marchetti, Soft Matter, 2023, 19, 7828–7835.
40. M. L. Blow, S. P. Thampi and J. M. Yeomans, Phys. Rev. Lett., 2014, 113, 248303.
41. H. Soni, W. Luo, R. A. Pelcovits and T. R. Powers, Soft Matter, 2019, 15, 6318–6330.
42. P. Gulati, F. Caballero, I. Kolvin, Z. You and M. C. Marchetti, Soft Matter, 2024, 20, 7703–7714.
43. R. Adkins, I. Kolvin, Z. H. You, S. Witthaus, M. C. Marchetti and Z. Dogic, Science, 2022, 377, 768–772.
44. R. Alert, C. Blanch-Mercader and J. Casademunt, Phys. Rev. Lett., 2019, 122, 088104.
45. L. Zhao, P. Gulati, F. Caballero, I. Kolvin, R. Adkins, M. C. Marchetti and Z. Dogic, Proc. Natl. Acad. Sci. U. S. A., 2024, 121, 2410345121.
46. S. Bhattacharyya and J. M. Yeomans, Phys. Rev. Lett., 2023, 130, 238201.
47. J. F. Joanny, F. Jülicher, K. Kruse and J. Prost, New J. Phys., 2007, 9, 422.
48. S. Bhattacharyya and J. M. Yeomans, Phys. Rev. E, 2024, 110, 024607.
49. A. Malevanets and J. M. Yeomans, Faraday Discuss., 1999, 112, 237–248.
50. A. Malevanets and J. M. Yeomans, Comput. Phys. Commun., 2000, 127, 105–112.
51. P. G. De Gennes and J. Prost, The Physics of Liquid Crystals, Oxford University Press, 1993.
52. A. N. Beris and B. J. Edwards, Thermodynamics of flowing systems: with internal microstructure, Oxford University Press, 1994.
53. S. Pokawanvit, Z. Chen, Z. You, L. Angheluta, M. C. Marchetti and M. J. Bowick, Phys. Rev. E, 2022, 106, 054610.
54. T. Krueger, H. Kusumaatmaja, A. Kuzmin, O. Shardt, G. Silva and E. Viggen, The Lattice Boltzmann Method: Principles and Practice, Springer, 2016.
55. R. C. V. Coelho, H. R. J. C. Figueiredo and M. M. Telo da Gama, Phys. Rev. Res., 2023, 5, 033165.
56. R. C. V. Coelho, N. A. M. Araújo and M. M. Telo da Gama, Philos. Trans. R. Soc., A, 2021, 379, 20200394.
57. S. Santhosh, M. R. Nejad, A. Doostmohammadi, J. M. Yeomans and S. P. Thampi, J. Stat. Phys., 2020, 180, 699–709.
58. F. G. Woodhouse and R. E. Goldstein, Phys. Rev. Lett., 2012, 109, 168105.
59. P. L. Townes and J. Holtfreter, J. Exp. Zool., 1955, 128, 53–120.
60. R. Moore, K. Q. Cai, D. O. Escudero and X.-X. Xu, Genesis, 2009, 47, 579–589.
61. O. Cochet-Escartin, T. T. Locke, W. H. Shi, R. E. Steele and E.-M. S. Collins, Biophys. J., 2017, 113, 2827–2841.
62. E. Batlle and D. G. Wilkinson, Cold Spring Harbor Perspect. Biol., 2012, 4, a008227.

---