2.1. Rate-of-deformation and vorticity tensors

A tensor of vital importance to fluid mechanics is the rate-of-deformation or strain-rate tensor $d_{\alpha \beta }$. The rate-of-deformation tensor is given by half of the time derivative of the metric on $\mathcal {M}_t$ (more formally, half of the Lie derivative with respect to the flow $\boldsymbol {V}$; Arroyo & DeSimone Reference Arroyo and DeSimone2009; Marsden & Hughes Reference Marsden and Hughes1994). If we take the time derivative of the induced metric, $g_{\alpha \beta }=\boldsymbol {e}_\alpha \boldsymbol {\cdot }\boldsymbol {e}_\beta$, as defined in (A1), then we find

(2.4)\begin{equation} \frac{\mathrm{d} g_{\alpha\beta}}{\mathrm{d}t} = \frac{\mathrm{d}\boldsymbol{e}_{\alpha}}{\mathrm{d}t}\boldsymbol{\cdot}\boldsymbol{e}_{\beta}+ \frac{\mathrm{d}\boldsymbol{e}_{\beta}}{\mathrm{d}t}\boldsymbol{\cdot}\boldsymbol{e}_{\alpha} = \boldsymbol{V}_{,\alpha}\boldsymbol{\cdot}\boldsymbol{e}_{\beta} + \boldsymbol{e}_{\alpha}\boldsymbol{\cdot}\boldsymbol{V}_{,\beta}. \end{equation}

Substituting for $\boldsymbol {V}_{,\alpha }=(v^\beta {}_{;\alpha } - v^{(n)} b^{\beta }{}_{\alpha })\boldsymbol {e}_{\beta } + (v^\beta b_{\alpha \beta } + v^{(n)}{}_{;\alpha })\hat {\boldsymbol {n}}$ (see (A13)) gives

(2.5)\begin{equation} d_{\alpha\beta}=\frac{1}{2}\,\frac{\mathrm{d} g_{\alpha\beta}}{\mathrm{d}t} = \frac{1}{2}\left(v_{\beta;\alpha} + v_{\alpha;\beta}\right) - v^{(n)} b_{\alpha\beta}. \end{equation}

Alternatively, one can start from the expression for the rate-of-deformation tensor in $\mathbb {R}^3$ and project back onto the manifold $\mathcal {M}_{t}$:

(2.6)\begin{equation} d_{\alpha\beta}\boldsymbol{e}^\alpha\otimes\boldsymbol{e}^\beta =\left[\tfrac{1}{2}\,\boldsymbol{V}_{,\alpha}\otimes\boldsymbol{e}^{\alpha} +\tfrac{1}{2}\left(\boldsymbol{V}_{,\alpha}\otimes\boldsymbol{e}^{\alpha}\right)^{\text{T}}\right]^{{\bullet}}, \end{equation}

where $({\cdot })^{\bullet }$ denotes projection onto $\mathcal {M}_t$ (i.e. contracting with tangent basis $\boldsymbol {e}_{\beta }$), which a quick calculation verifies gives the same result. Note that $\boldsymbol {V}\in \mathbb {R}^3$ but $\boldsymbol {V}$ has support only on $\mathcal {M}_t$, and that the formal derivation of this requires some subtle constructions such as extensions of vector fields (Lee Reference Lee1997).

The vorticity tensor can be derived in an equivalent way, giving a result identical to that from the flat geometry version:

(2.7)\begin{equation} \varOmega_{\alpha\beta}=\tfrac{1}{2}\left(v_{\beta;\alpha}-v_{\alpha;\beta}\right), \end{equation}

which can be seen by noting that the terms of the antisymmetrised deformation tensor proportional to $b_{\alpha \beta }$ vanish by symmetry under exchange of indices of the second fundamental form.

Here, we will consider incompressible fluids, the condition for which is that the rate-of-deformation tensor is traceless:

(2.8)\begin{equation} d^\alpha{}_\alpha = \boldsymbol{\nabla}_\alpha v^\alpha - 2 H v^{(n)}=0. \end{equation}

In the absence of normal flows, this reproduces the classic incompressibility condition $\boldsymbol {\nabla }_\alpha v^\alpha =0$.

2.2. Objective rate of a vector field

The subject of objective rates is contentious in continuum mechanics, especially where changing geometry is concerned (Marsden & Hughes Reference Marsden and Hughes1994; Nitschke & Voigt Reference Nitschke and Voigt2022). The challenge, in this case, is to ensure that tangential flows are Eulerian – i.e. coordinates are not advected with the flow – yet also properly account for out-of-plane surface dynamics, which are necessarily Lagrangian – i.e. coordinates are advected with surface movement (Waxman Reference Waxman1984). We adopt a pragmatic approach, and choose a projected co-rotational rate. That is, we take the co-rotational derivative in $\mathbb {R}^3$, and project down onto the surface so as to recover the co-rotational objective rates for fixed-curved and flat surfaces, while still accounting for Lagrangian deformations of the surface in the normal direction.

For a vector field $\boldsymbol {T}$, which is confined to the tangent bundle of $\mathcal {M}_{t}$ – i.e. $\boldsymbol {T}=T^\alpha \boldsymbol {e}_\alpha \in \mathcal {T}(\mathcal {M}_t)$ – this takes the form

(2.9)\begin{equation} {D}_t\boldsymbol{T}^{{\bullet}} =\left[\frac{\mathrm{d}\boldsymbol{T}}{\mathrm{d}t} + \varOmega\boldsymbol{\cdot}\boldsymbol{T}\right]^{{\bullet}} = \left[\partial_t\boldsymbol{T} +v^\alpha(\boldsymbol{T})_{,\alpha} + \varOmega\boldsymbol{\cdot}\boldsymbol{T}\right]^{{\bullet}}. \end{equation}

The first term is given by

(2.10)\begin{align} \partial_t\boldsymbol{T} &= \partial_{t}\left(T^\alpha\boldsymbol{e}_\alpha\right) = \partial_tT^\alpha \boldsymbol{e}_\alpha + T^\alpha \partial_t\boldsymbol{e}_\alpha\nonumber\\ & = \left(\partial_tT^\alpha - v^{(n)} T^\beta b_{\beta}{}^{\alpha}\right)\boldsymbol{e}_\alpha + T^\alpha v^{(n)}{}_{;\alpha} \hat{\boldsymbol{n}}, \end{align}

where we have included only contributions from changing the basis vectors due to normal velocities, since we are in an Eulerian frame tangentially and therefore tangential velocities do not advect coordinates. The second term is given by

(2.11)\begin{equation} v^{\alpha}\,\boldsymbol{\nabla}_\alpha(\boldsymbol{T}) = v^\alpha\left(T^\beta \boldsymbol{e}_\beta\right)_{,\alpha} = v^\alpha T^\beta{}_{;\alpha}\boldsymbol{e}_\beta + T^\beta v^{\alpha} b_{\alpha\beta}\hat{\boldsymbol{n}}, \end{equation}

and finally the third term is simply given by

(2.12)\begin{equation} \varOmega\boldsymbol{\cdot}\boldsymbol{T}=\varOmega^{\alpha}{}_{\beta}T^\beta\boldsymbol{e}_{\alpha}. \end{equation}

Putting this all together and projecting back onto the manifold gives

(2.13)\begin{equation} {D}_tT^\alpha = \partial_tT^\alpha + v^\beta T^\alpha{}_{;\beta} + \varOmega^{\alpha}{}_{\beta}T^\beta - v^{(n)} T^\beta b_{\beta}{}^{\alpha}, \end{equation}

which is identical to the equivalent expression in flat geometry, with the exception of the final term.

To understand the origin of this term, consider a tube in cylindrical polar coordinates that is expanding outwards with velocity $v^{(n)}=\dot {R}$, where $R$ is the radius. There is no tangential flow, and the director field is assumed to be of constant magnitude in the $\boldsymbol {e}_{\theta }$ direction (i.e. $\boldsymbol {T}=T\boldsymbol {e}_{\theta }$). In this case, the second fundamental form has one non-zero component in the $\boldsymbol {e}_\theta$ direction. As a result, as $R$ increases, so does the size of the basis vector $\boldsymbol {e}_\theta$, and we require that ${D}_tT^{\theta } = \partial _t T + \dot {R} T/R=0$ for the objective rate to be zero (fixed size for $\boldsymbol {T}$), which gives $T=T_0\exp (-\dot {R}/R)$.

2.3. Elastic forces

In general, we will consider surfaces with free energies of the form

(2.14)\begin{equation} \mathcal{F}=\int_{\mathcal{M}_t}F\left(g_{\alpha\beta},b_{\alpha\beta},T_{\alpha}\right)\mathrm{d}A, \end{equation}

where $\mathrm {d}A=\sqrt {|g|}\,\mathrm {d}\kern0.7pt x^1\,\mathrm {d}\kern0.7pt x^2$ is the area element on the manifold. Whilst other free energies fit into this class, for example the tilt field of a lipid monolayer (Selinger, MacKintosh & Schnur Reference Selinger, MacKintosh and Schnur1996; Terzi & Deserno Reference Terzi and Deserno2017), we will focus on the concrete example of the free energy of a nematic liquid crystal. Here, we will make use of standard descriptions of liquid crystals to define the molecular field as

(2.15)\begin{equation} h_\alpha ={-} \frac{\delta \mathcal{F}}{\delta T^\alpha}, \end{equation}

which we will make use of in order to derive the hydrodynamics of the ordered surface. The force given by varying the surface shape is then given by

(2.16)\begin{equation} \boldsymbol{f}_{{e}}=f^\alpha_{{e}}\boldsymbol{e}_\alpha + f^{(n)}_{{e}}\hat{\boldsymbol{n}} ={-}\frac{\delta \mathcal{F}}{\delta \boldsymbol{R}}. \end{equation}

In the case of nematic liquid crystals, the one-constant approximation to the Frank free energy (Frank Reference Frank1958; DeGennes & Prost Reference DeGennes and Prost1993; Napoli & Vergori Reference Napoli and Vergori2012, Reference Napoli and Vergori2016) is given by

(2.17)\begin{equation} \mathcal{F}^{{LC}} = \int_{\mathcal{M}_t} \frac{\mathcal{K}}{2}\,|\tilde{\boldsymbol{\nabla}}_{s}\boldsymbol{T}|^2\,\mathrm{d}A , \end{equation}

where $\boldsymbol {T} = T^\beta \boldsymbol {e}_\beta$ is the nematic order parameter (we choose to use $T^\beta$ here rather than the usual $n^\beta$ in the nematic liquid crystal literature, in order to avoid confusion with the surface normal vector $\hat {\boldsymbol {n}}$) on the tangent bundle $\mathcal {T}(\mathcal {M}_t)$, and $\tilde {\boldsymbol {\nabla }}_s$ is the surface derivative (not the covariant derivative) (Napoli & Vergori Reference Napoli and Vergori2012). The reason for our choice of the surface derivative, as opposed to the covariant derivative, is so as to account correctly for the bend due to extrinsic curvature. This is likely an important coupling for true active nematic surfaces, but is perhaps less relevant for nematic-like material surfaces such as epithelial tissues. The surface derivative can be written in terms of the covariant derivative and second fundamental form as

(2.18)\begin{equation} \tilde{\boldsymbol{\nabla}}_s \boldsymbol{T} = \boldsymbol{\nabla}_\alpha T^\beta \boldsymbol{e}^\alpha\otimes\boldsymbol{e}_\beta + T^\beta b_{\alpha\beta} \boldsymbol{e}^\alpha\otimes\hat{\boldsymbol{n}}, \end{equation}

thus the full free energy is given by

(2.19)\begin{equation} \mathcal{F}^{{LC}} = \int_{\mathcal{M}_t} \frac{\mathcal{K}}{2}\left[\boldsymbol{\nabla}_\alpha T^\beta\,\boldsymbol{\nabla}^\alpha T_\beta + T^\beta T_\gamma b_{\alpha\beta}b^{\alpha\gamma}\right]\mathrm{d}A. \end{equation}

We will now compute the associated functional derivatives in the free energy (${\mathcal {F}}$). The functional derivative with respect to $T^\alpha$ gives the negative of the molecular field $h_\alpha$ as

(2.20)\begin{equation} \frac{\delta \mathcal{F}^{{LC}}}{\delta T^\alpha} ={-} h_\alpha ={-}\mathcal{K} \left[ \Delta T_\alpha + \left(K\delta_\alpha^\gamma - 2H b_\alpha{}^\gamma\right)T_\gamma \right]. \end{equation}

We will make use of the rate formulas calculated in Appendix A to perform functional variations of the form

(2.21)\begin{equation} \left. \begin{aligned} \delta g_{\alpha\beta} = \delta t\,\frac{\mathrm{d}}{\mathrm{d}t}\left(g_{\alpha\beta}\right)+O(\delta t^2),\cr \delta b_{\alpha\beta} = \delta t\,\frac{\mathrm{d}}{\mathrm{d}t}\left(b_{\alpha\beta}\right)+O(\delta t^2),\cr \delta \varGamma^{\gamma}_{\alpha\beta} = \delta t\,\frac{\mathrm{d}}{\mathrm{d}t}\left(\varGamma^\gamma_{\alpha\beta}\right)+O(\delta t^2). \end{aligned} \right\} \end{equation}

Using this, we find

(2.22) $$\begin{align} \delta\mathcal{F} &= \int_{\mathcal{M}_t} \frac{\mathcal{K}}{2} \left[\vphantom{\left( T^\beta{}_{;\alpha} T_\beta{}^{;\alpha} + T^\beta T_\gamma b_{\alpha\beta} b^{\alpha\gamma}\right)} 2 \delta(\varGamma^\beta_{\alpha\mu}) T^\mu g^{\alpha\gamma} g_{\beta\delta} T^{\delta}{}_{;\gamma} + T^\beta{}_{;\alpha} \delta(g^{\alpha\gamma}) g_{\beta\delta} T^{\delta}{}_{;\gamma} + T^\beta{}_{;\alpha} g^{\alpha\gamma} \delta(g_{\beta\delta}) T^{\delta}{}_{;\gamma}\right.\nonumber\\ &\quad\left. {}+ 2 \delta (b_{\alpha\beta}) g^{\alpha\delta}b_{\delta\gamma} T^\gamma T^\beta + T^\beta T^\gamma b_{\alpha\beta}b_{\delta\gamma} \delta(g^{\alpha\delta}) -2H v^{(n)} \delta t\left( T^\beta{}_{;\alpha} T_\beta{}^{;\alpha} + T^\beta T_\gamma b_{\alpha\beta} b^{\alpha\gamma}\right)\right]\mathrm{d}A. \end{align}$$

Focusing purely on the normal variation, we find

(2.23) \begin{align} \delta \mathcal{F}&= \int_{\mathcal{M}_t}\mathcal{K} \delta t\left\{ \left[ (v^{(n)} b_{\alpha}{}^\gamma)_{;\beta} - (v^{(n)}b_{\alpha\beta})^{;\gamma} - (v^{(n)}b^\gamma{}_\beta)_{;\alpha}\right]T^\alpha T^\beta{}_{;\gamma} + T^\beta{}_{;\alpha}T_{\beta;\gamma}b^{\alpha\gamma}v^{(n)} \right.\nonumber\\ &\quad - T^\beta{}_{;\alpha} T^{\gamma;\alpha}b_{\gamma\beta}v^{(n)} +(b^\alpha{}_\gamma T^\gamma T^\beta)_{;\beta\alpha}v^{(n)} - T^\gamma T^\beta \left(2H b_{\alpha\beta}b^\alpha{}_\gamma - K b_{\beta\gamma} \right)v^{(n)}\nonumber\\ &\left.\quad+\,T^\beta T^\gamma b_{\alpha\beta} b_{\delta\gamma}b^{\alpha\delta}v^{(n)}- H \left[ T^\beta{}_{;\alpha} T_\beta{}^{;\alpha} + T^\beta T_\gamma b_{\alpha\beta} b^{\alpha\gamma}\right]v^{(n)}\right\}\mathrm{d}A. \end{align}

Notice that the first terms from the variation of the covariant derivative include derivatives in the normal velocity. Integrating these by parts and pushing terms to the boundary, we can find the following functional variation in the bulk of the surface:

(2.24)\begin{align} \delta\mathcal{F}& = \int_{\mathcal{M}_t}\mathcal{K} v^{(n)}\delta t\left\{ b_{\alpha\beta}(T^\alpha T^\beta{}_{;\gamma})^{;\gamma} + b^\gamma{}_\beta(T^\alpha T^\beta{}_{;\gamma})_{;\alpha}- b_{\alpha}{}^\gamma (T^\alpha T^\beta{}_{;\gamma})_{;\beta} + T^\beta{}_{;\alpha}T_{\beta;\gamma}b^{\alpha\gamma}\right.\nonumber\\ &\quad - T^\beta{}_{;\alpha} T^{\gamma;\alpha}b_{\gamma\beta} + (b^\alpha{}_\gamma T^\gamma T^\beta)_{;\beta\alpha} - T^\gamma T^\beta \left(2H b_{\alpha\beta}b^\alpha{}_\gamma - K b_{\beta\gamma} \right) + T^\beta T^\gamma b_{\alpha\beta} b_{\delta\gamma}b^{\alpha\delta}\nonumber\\ &\left.\quad -\, H \left[ T^\beta{}_{;\alpha} T_\beta{}^{;\alpha} + T^\beta T_\gamma b_{\alpha\beta} b^{\alpha\gamma}\right]\right\}\,\mathrm{d}A + \text{boundary terms}. \end{align}

From this, and noting that the second and third terms cancel due to symmetry of $T^\alpha T^\beta$ and $b_{\alpha \beta }$, we find the forces in the normal direction from the free energy as

(2.25)\begin{align} &f^{(n)}_{{e}}={-}\mathcal{K}\left\{ b_{\alpha\beta}(T^\alpha T^\beta{}_{;\gamma})^{;\gamma} + T^\beta{}_{;\alpha}\left(T_{\beta;\gamma}b^{\alpha\gamma} - T^{\gamma;\alpha}b_{\gamma\beta}\right)+(b^\alpha{}_\gamma T^\gamma T^\beta)_{;\beta\alpha}\right.\nonumber\\ &\left.\quad {}- T^\gamma T^\beta \left(2H b_{\alpha\beta}b^\alpha{}_\gamma - K b_{\beta\gamma} \right) + T^\beta T^\gamma b_{\alpha\beta} b_{\delta\gamma}b^{\alpha\delta} - H \left[ T^\beta{}_{;\alpha} T_\beta{}^{;\alpha} + T^\beta T_\gamma b_{\alpha\beta} b^{\alpha\gamma}\right]\right\}. \end{align}

Note that the third term gives forces that go like the second derivative of the curvature tensor, but traced out along the directions of the nematic order parameter. Such a term is essentially an anisotropic bending rigidity that would not be found in the case where the covariant derivative, rather than surface derivative, was used in the free energy (Napoli & Vergori Reference Napoli and Vergori2012; Santiago, Chacón-Acosta & Monroy Reference Santiago, Chacón-Acosta and Monroy2019; Santiago & Monroy Reference Santiago and Monroy2020). This gives equations that are significantly different to many of those currently analysed in the literature, both equilibrium and dynamical (Frank & Kardar Reference Frank and Kardar2008; Santiago Reference Santiago2018; Hoffmann et al. Reference Hoffmann, Carenza, Eckert and Giomi2022; Salbreux et al. Reference Salbreux, Jülicher, Prost and Callen-Jones2022).

Similarly for the tangential variations, we find the elastic force

(2.26)$$\begin{gather} f^{\alpha}_{{e}} ={-}\frac{\mathcal{K}}{2} \left\{ \left(T^\beta T^{\alpha;\gamma}\right)_{;\beta\gamma} + \left(T_{\beta;\gamma} T^{\beta;\alpha}\right)^{;\gamma} - \left(T_{\beta;\gamma} T^{\alpha;\gamma}\right)^{;\beta} - 2 \left(b^\alpha{}_\beta b^\gamma{}_\mu T^{\mu\beta}\right)_{;\gamma} \right.\nonumber\\ \left.{}- 2 \left( b^\gamma{}_\mu T^{\mu\beta}\right)_{;\beta} b^\alpha{}_\gamma + \left(T^\beta T^\gamma b_{\mu\beta} b^\alpha{}_{\gamma}\right)^{;\mu}\right\}. \end{gather}$$

This is essentially a generalised version of the Ericksen force (DeGennes & Prost Reference DeGennes and Prost1993). For simplicity, we will neglect this tangential component in our later examples as it leaves the phenomenology of the equations unchanged. This approximation is taken in many exact calculations in active hydrodynamics (Hoffmann et al. Reference Hoffmann, Carenza, Eckert and Giomi2022; Giomi et al. Reference Giomi, Bowick, Mishra, Sknepnek and Marchetti2014; Kruse et al. Reference Kruse, Joanny, Jülicher, Prost and Sekimoto2004; Edwards & Yeomans Reference Edwards and Yeomans2009); however, in principle this term should be included in full nonlinear solutions.

2.4. Surface tension and director constraints

We consider two constraints upon our system whose contributions to the dynamical equations can be derived generically from a Rayleigh dissipation functional (DeGennes & Prost Reference DeGennes and Prost1993; Doi Reference Doi2011). However, we avoid a full variational treatment here for the sake of brevity, and instead simply state the contributions whose derivation can be found elsewhere (DeGennes & Prost Reference DeGennes and Prost1993; Arroyo & DeSimone Reference Arroyo and DeSimone2009; Marchetti et al. Reference Marchetti, Joanny, Ramaswamy, Liverpool, Prost, Rao and Simha2013).

The first constraint is that of constant area, or incompressibility, which is expressed by (2.8). The corresponding Lagrange multiplier is the surface tension $\zeta$, which acts isotropically via stresses $\sigma ^{(st)}_{\alpha \beta }=\zeta g_{\alpha \beta }$. Taking the surface divergence of $\sigma ^{(st)}$ gives the following contribution to the force balance condition:

(2.27)\begin{equation} \boldsymbol{f}^{c} = \boldsymbol{\nabla}^\alpha \zeta\boldsymbol{e}_\alpha + 2H \zeta \hat{\boldsymbol{n}}. \end{equation}

The second constraint is to fix the magnitude of the nematic director (DeGennes & Prost Reference DeGennes and Prost1993). This is achieved by imposing

(2.28)\begin{equation} T_\alpha {D}_t T^\alpha = 0, \end{equation}

which implies $T^\alpha T_\alpha =\textrm {constant}$, with the constant specified by the initial conditions on $T^\alpha$ (Kruse et al. Reference Kruse, Joanny, Jülicher, Prost and Sekimoto2005). We label the associated Lagrange multiplier by $\lambda$, and the corresponding contribution to the director current (Prost, Jülicher & Joanny Reference Prost, Jülicher and Joanny2015) is just

(2.29)\begin{equation} J^\alpha_{c} ={-}\lambda T^\alpha. \end{equation}

This term can be viewed from a variational perspective as either a constraint on the Rayleigh dissipation functional of the form $\int \lambda T_\alpha {D}_t T^\alpha \,\mathrm {d}A$, or a contribution to the molecular field from a constraint on the free energy of the form $\int \frac {1}{2}\lambda \gamma T_\alpha T^\alpha \,\mathrm {d}A$ (DeGennes & Prost Reference DeGennes and Prost1993).

2.5. Constitutive relation and dissipative forces

The viscous stress tensor in the tangent bundle is given as an Onsager expansion (Doi Reference Doi2011; DeGennes & Prost Reference DeGennes and Prost1993) in the thermodynamic fluxes. At lowest non-trivial order, it is given by

(2.30)\begin{equation} \sigma_{\alpha\beta} = 2 \eta d_{\alpha\beta} + \frac{\nu}{2}\left(T_\alpha h_\beta + T_\beta h_\alpha\right), \end{equation}

where $h_\alpha$ is the molecular field (2.20), $d_{\alpha \beta }$ is the rate-of-deformation tensor, $\eta$ is the shear viscosity, and $\nu$ gives the spin connection coefficients of the vector field (essentially the ratio of anisotropic viscosity to isotropic viscosity DeGennes & Prost Reference DeGennes and Prost1993). Here, stable solutions require $|\nu |>1$, with $\nu <-1$ corresponding to rod-like nematics, and $\nu >1$ to disk-like nematics (Edwards & Yeomans Reference Edwards and Yeomans2009; Marchetti et al. Reference Marchetti, Joanny, Ramaswamy, Liverpool, Prost, Rao and Simha2013). Note that as we consider only incompressible fluids, we discard the bulk viscosity terms and associated spin connections (DeGennes & Prost Reference DeGennes and Prost1993). The dissipative forces on our surface are then given as

(2.31)\begin{equation} \boldsymbol{f}^{d} = \boldsymbol{\nabla}_\beta \sigma^{\alpha\beta}\boldsymbol{e}_\alpha + \sigma^{\alpha\beta}b_{\alpha\beta} \hat{\boldsymbol{n}}. \end{equation}

For the terms in the tangent bundle we find

(2.32)\begin{equation} \boldsymbol{\nabla}_\beta\sigma^{\alpha\beta} = 2 \eta d^{\alpha\beta}{}_{;\beta} + \frac{\nu}{2}\left(T^\alpha h^\beta + T^\beta h^\alpha\right)_{;\beta}, \end{equation}

where we can unpack the first term to

(2.33)\begin{equation} 2\eta d^{\alpha\beta}{}_{;\beta} = \eta \left(v^{\alpha;\beta} + v^{\beta;\alpha} - 2 v^{(n)}b^{\alpha\beta}\right)_{;\beta} = \eta\varDelta v^\alpha + \eta v^\beta{}_{;\alpha\beta} - 2\eta \left(v^{(n)}b^{\alpha\beta}\right)_{;\beta}, \end{equation}

where $\varDelta =\boldsymbol {\nabla }_\alpha \boldsymbol {\nabla }^\alpha$ is the Bochner Laplacian. By making use of the Codazzi equation (A9) and the contracted Gauss equation (A12), we find

(2.34)\begin{equation} 2\eta d^{\alpha\beta}{}_{;\beta} = \eta\left[ \varDelta v^\alpha - 2\left(v^{(n)}H\right)^{;\alpha} + K v^{\alpha} + 2\left(2H g^{\alpha\beta} - b^{\alpha\beta}\right)v^{(n)}{}_{;\beta}\right]. \end{equation}

Thus the tangential forces are given by

(2.35)\begin{align} \sigma^{\alpha\beta}{}_{;\beta} &= \frac{\nu}{2}\left(T^\alpha h^\beta + T^\beta h^\alpha\right)_{;\beta}\nonumber\\ &\quad +\eta\left[ \varDelta v^\alpha - 2\left(v^{(n)}H\right)^{;\alpha} + K v^{\alpha} + 2\left(2H g^{\alpha\beta} - b^{\alpha\beta}\right)v^{(n)}{}_{;\beta}\right] . \end{align}

The normal forces are given by contracting the stress tensor with the second fundamental form

(2.36)\begin{equation} \sigma^{\alpha\beta}b_{\alpha\beta} = 2\eta b^{\alpha\beta}v_{\beta;\alpha} - 4\eta(2 H^2 - K)v^{(n)} + \nu T^{\alpha}h^{\beta}b_{\alpha\beta}. \end{equation}

2.6. Broken symmetry terms

We can also include broken symmetry terms in the stress tensor that contribute nothing to the overall dissipation as their dynamics are associated with Goldstone modes in the system, hence the name ‘broken symmetry’ (DeGennes & Prost Reference DeGennes and Prost1993; Kruse et al. Reference Kruse, Joanny, Jülicher, Prost and Sekimoto2005). For vector fields, these contribute an antisymmetric part to the stress tensor of the form

(2.37)\begin{equation} \sigma^{(bs)}_{\alpha\beta} = \tfrac{1}{2}\left(T_\alpha h_\beta - T_\beta h_\alpha\right). \end{equation}

Note that these terms contribute to only the tangential forces as they vanish in the normal due to the symmetry of the second fundamental form upon index relabelling. The broken symmetry forces are thus given by

(2.38)\begin{equation} \boldsymbol{f}^{bs} = \boldsymbol{\nabla}_\beta \sigma^{(bs)\alpha\beta}\boldsymbol{e}_\alpha + \sigma^{(bs)\alpha\beta}b_{\alpha\beta} \hat{\boldsymbol{n}} = \tfrac{1}{2}\left(T^\alpha h^\beta - T^\beta h^\alpha\right)_{;\beta}\boldsymbol{e}_\alpha. \end{equation}

2.7. Director terms

Expanding out in an Onsager fashion for the director rate in terms of nematic relaxation, we find the couplings

(2.39)\begin{equation} J^{\alpha}_{T} = \gamma^{{-}1} h^\alpha -\nu d^\alpha{}_{\beta}T^\beta, \end{equation}

which give terms in a similar manner to the flat case (Kruse et al. Reference Kruse, Joanny, Jülicher, Prost and Sekimoto2005). Note that the constant $\nu$ is the same as for the director terms in the stress tensor by the Onsager reciprocal relations. Adding this together with the constraint current (2.29) gives the full Onsager expansion for the director rate as

(2.40)\begin{equation} {D}_t T^\alpha =\gamma^{{-}1} h^\alpha - \nu d^\alpha{}_\beta T^\beta - \lambda T^\alpha. \end{equation}

2.8. Active terms

We consider a simple extension to the tangential stress tensor that represents the dipolar-like active forces that are characteristic of many active nematic systems (Sanchez et al. Reference Sanchez, Chen, DeCamp, Heymann and Dogic2012; Marchetti et al. Reference Marchetti, Joanny, Ramaswamy, Liverpool, Prost, Rao and Simha2013). For a detailed introduction to active matter, we refer the reader to the reviews Marchetti et al. (Reference Marchetti, Joanny, Ramaswamy, Liverpool, Prost, Rao and Simha2013) and Ramaswamy (Reference Ramaswamy2017), and for active liquid crystals in particular to Doostmohammadi et al. (Reference Doostmohammadi, Ignés-Mullol, Yeomans and Sagués2018).

Active liquid crystals have generated significant interest in recent years due to several experimental realisations (Sanchez et al. Reference Sanchez, Chen, DeCamp, Heymann and Dogic2012; Keber et al. Reference Keber, Loiseau, Sanchez, DeCamp, Giomi, Bowick, Marchetti, Dogic and Bausch2014) that have revealed the fascinating interplay between active stresses and geometry of the director, particularly in the coupling of geometry of the surfaces to the dynamics of topological defects (Giomi et al. Reference Giomi, Bowick, Mishra, Sknepnek and Marchetti2014; Khoromskaia & Alexander Reference Khoromskaia and Alexander2017; Pearce et al. Reference Pearce, Ellis, Fernandez-Nieves and Giomi2019; Pearce Reference Pearce2020). In addition such systems have been suggested as minimal models for tissues (Blanch-Mercader et al. Reference Blanch-Mercader, Guillamat, Roux and Kruse2021a,Reference Blanch-Mercader, Guillamat, Roux and Kruseb) and the actomyosin cell cortex (Kruse et al. Reference Kruse, Joanny, Jülicher, Prost and Sekimoto2005).

Here, we consider the standard form for the active nematic stress tensor:

(2.41)\begin{equation} \sigma^{{(a)}}_{\alpha\beta} = \varSigma Q_{\alpha\beta} = \varSigma \left(T_\alpha T_\beta - \tfrac{1}{2} g_{\alpha\beta}\right). \end{equation}

For $\varSigma >0$, such stresses are contractile in nature, and for $\varSigma <0$, they are extensile. This form represents the simplest form of anisotropic active stress; however, more complicated couplings are possible (Naganathan et al. Reference Naganathan, Fürthauer, Nishikawa, Jülicher and Grill2014; Salbreux & Jülicher Reference Salbreux and Jülicher2017; Salbreux et al. Reference Salbreux, Jülicher, Prost and Callen-Jones2022). The active forces from (2.41) are given by

(2.42)\begin{equation} \boldsymbol{f}^a = \sigma^{{(a)}\alpha\beta}{}_{;\beta} \boldsymbol{e}_{\alpha}+ \sigma^{{(a)}\alpha\beta}b_{\alpha\beta}\hat{\boldsymbol{n}}. \end{equation}

2.9. Full dynamical equations

Putting together all our force balance terms, constraints and polarisation currents gives the following system of equations, which are essentially ordered fluid versions of the equations derived for isotropic fluid membranes in Arroyo & DeSimone (Reference Arroyo and DeSimone2009).

(2.43)\begin{align} &\text{Polarisation dynamics:}\nonumber\\ &\quad {D}_t T^\alpha =\gamma^{{-}1} h^\alpha - \nu d^\alpha{}_\beta T^\beta - \lambda T^\alpha.\end{align}

(2.44)\begin{align} &\text{Tangential force balance:}\nonumber\\ &\quad \frac{1}{2}\left(T^\alpha h^\beta - T^\beta h^\alpha\right)_{;\beta} +\eta\left[ \varDelta v^\alpha - 2\left(v^{(n)}H\right)^{;\alpha} + K v^{\alpha} + 2\left(2H g^{\alpha\beta} - b^{\alpha\beta}\right)v^{(n)}{}_{;\beta}\right]\nonumber\\ &\qquad + \frac{\nu}{2}\left(T^\alpha h^\beta + T^\beta h^\alpha\right)_{;\beta} + f^\alpha_{{e}} + \zeta^{;\alpha} + \sigma^{{(a)}\alpha\beta}{}_{;\beta}=0. \end{align}

(2.45)\begin{align} &\text{Normal force balance:}\nonumber\\ &\quad 2\eta b^{\alpha\beta}v_{\beta;\alpha} - 4\eta(2 H^2 - K )v^{(n)} + \nu T^{\alpha}h^{\beta}b_{\alpha\beta} + 2 H\zeta + f^{(n)}_{{e}} + \sigma^{{(a)}\alpha\beta}b_{\alpha\beta}=0. \end{align}

(2.46)\begin{align} &\text{Constraint equations (constant director and incompressibility):}\nonumber\\ &\quad T_\alpha{D}_t T^\alpha=0,\quad v^\alpha{}_{;\alpha}-2H v^{(n)}=0.\end{align}

(2.47)\begin{align} &\text{Closure between surface velocity and geometry:}\nonumber\\ &\quad \boldsymbol{V}= \frac{\mathrm{d}\boldsymbol{R}}{\mathrm{d}t}. \end{align}

2.10. Dimensionless numbers

There is an interesting point to make regarding the normal force balance equation. It contains the term $\eta b^{\alpha \beta }v_{\beta ;\alpha }$, which describes the curvature-induced coupling between viscous forces and shape. Such couplings have been discussed extensively for isotropic fluid membranes (Al-Izzi, Sens & Turner Reference Al-Izzi, Sens and Turner2020a; Sahu et al. Reference Sahu, Glisman, Tchoufag and Mandadapu2020a). These phenomena can be understood in terms of the dimensionless Scriven–Love number (Scriven Reference Scriven1960)

(2.48)\begin{equation} {SL} = \frac{O(\eta b^{\alpha\beta}v_{\beta;\alpha})}{O(\text{bending forces})} = \frac{\eta V L}{\mathcal{K}}, \end{equation}

where $V$ and $L$ are characteristic velocities and length scales, and we have assumed that $\textrm {bending forces}\sim \mathcal {K} L^{-3}$ (i.e., like the bending energy and Laplacian of the mean curvature).

In addition, however, there is a new anisotropic coupling due to terms $\sim \nu T^{\alpha }h^{\beta }b_{\alpha \beta }$ which, when compared to bending forces, could be thought of as a ‘liquid crystalline Scriven–Love number’:

(2.49)\begin{equation} {SL}^{(LC)} = \frac{O(\nu T_{\alpha}h_{\beta} b^{\alpha\beta})}{O(\text{bending forces})} = \nu. \end{equation}

As such, the dimensionless number associated with out-of-plane liquid crystalline dissipative forces is just the spin connection coefficient $\nu$.

In addition to this, there is also the dimensionless active stress. As this is a number comparing tangential stress to bending stress, we can consider it as an analogue of the Föppl–von Kármán number ${Y}$ that would typically compare stretching to bending forces in passive shell theory (Audoly & Pomeau Reference Audoly and Pomeau2010). In our case, this is given by

(2.50)\begin{equation} {Y}^{({Active})} = \frac{\varSigma L^2}{\mathcal{K}}. \end{equation}

Note that this implies that in highly curved geometries, bending forces will suppress active deformations.

3.1. Relaxation dynamics about a perturbed tube

We start with a simple example of a tube of radius $r(\theta,z,t)= r + \epsilon \,u(\theta,z,t)$, where $\epsilon$ is some dimensionless parameter that will be considered small. Such tubes have been shown in the case of isotropic membranes to have non-negligible Scriven–Love effects (Sahu et al. Reference Sahu, Glisman, Tchoufag and Mandadapu2020a). The surface is prescribed by the vector $\boldsymbol {R}=(r\cos \theta,r\sin \theta,z)$. We can thus calculate the induced basis and metric

(3.1)\begin{equation} [g_{\alpha\beta}]=\left(\begin{matrix} r(r+2\epsilon u) & 0\\ 0 & 1 \end{matrix}\right)+O(\epsilon^2), \end{equation}

and the second fundamental form

(3.2)\begin{equation} [b_{\alpha\beta}]=\left(\begin{matrix} \epsilon u_{,\theta\theta}-\epsilon u-r & \epsilon u_{,\theta z }\\ \epsilon u_{,\theta z} & \epsilon u_{,zz} \end{matrix}\right)+O(\epsilon^2). \end{equation}

This gives mean and Gaussian curvature to first order as

(3.3)\begin{gather} H=\frac{\epsilon}{2}\left( u_{,zz}+\frac{1}{r^2}\,u_{,\theta\theta}+ \frac{1}{r^2}\,u\right)-\frac{1}{2 r}, \end{gather}

(3.4)\begin{gather} K={-}\frac{\epsilon u_{,zz}}{r}. \end{gather}

We write the director field as $\boldsymbol {T}=({\epsilon T^\theta }/{r})\boldsymbol {e}_\theta + (1+\epsilon T^z) \boldsymbol {e}_z +O(\epsilon ^2)$, where we have chosen the $O(1)$ term to be the unit vector in the $z$ direction such that the molecular field vanishes for the unperturbed state (Napoli & Vergori Reference Napoli and Vergori2012). Thus at lowest order, the constraint equation on director unity becomes

(3.5)\begin{equation} \partial_t T^z=0, \end{equation}

and all the dynamics are in the $\theta$ direction, i.e. perpendicular to the initial zeroth order nematic ordering, as one would expect for the perturbative dynamics of a unit vector field. We write the velocity as $\boldsymbol{V} =\epsilon (v^\theta /r\boldsymbol {e}_\theta +v^z\boldsymbol {e}_z + \partial _t u \hat {\boldsymbol {n}})$, so the continuity equation becomes

(3.6)\begin{equation} \frac{\partial_t u}{r}+\frac{v^\theta_{,\theta}}{r}+v^z_{,z}=0. \end{equation}

The molecular field is given by

(3.7)\begin{equation} h^\alpha\boldsymbol{e}_\alpha = \mathcal{K} \epsilon \left(T^\theta_{,zz} + \frac{1}{r^2}\,T^\theta_{,\theta\theta} - \frac{1}{r^2}\,T^\theta + \frac{2}{r^2}\,u_{,\theta z}\right)\frac{1}{r}\,\boldsymbol{e}_\theta + \mathcal{K} \epsilon \left(\frac{T^z_{,\theta\theta}}{r^2} + T^z_{,zz}\right)\boldsymbol{e}_z, \end{equation}

and can be substituted to then find the polarisation rate equation whose components are given by

(3.8)\begin{align} &\boldsymbol{e}_\theta:\quad \epsilon \left[- \frac{\mathcal{K} }{\gamma r} \left(T^\theta_{,zz} + \frac{T^\theta_{,\theta\theta}}{r^2} - \frac{T^\theta}{r^2} + 2\,\frac{u_{,\theta z}}{r^2}\right)+\frac{1}{r} \left( \partial_t T^\theta + \frac{1}{2}\,v^\theta_{,z} +\frac{1}{2r}\, v^z_{,\theta}\right)\right.\nonumber\\ &\left.\qquad\quad{}+ \frac{1}{r} \left( \lambda T^\theta + \frac{\nu}{2}\,v^\theta_{,z} + \frac{\nu}{2r}\,v^z_{,\theta}\right)\right]=0, \end{align}

(3.9)\begin{align} &\boldsymbol{e}_z:\quad \lambda-\epsilon\left[\frac{\mathcal{K} }{\gamma} \left(\frac{1}{r^2}\,T^z_{,\theta\theta}+ T^z_{,zz}\right) + \partial_t T^z + \lambda T^z + \nu v^z_{,z}\right]=0. \end{align}

The modified covariant Stokes equation (tangential force balance) is given componentwise by

(3.10)\begin{align} &\boldsymbol{e}_\theta:\quad \frac{\zeta_{,\theta} (r-2 \epsilon u)}{r^3}+\frac{\mathcal{K} (\nu -1) \epsilon}{2 r^3} \left(r^2 T^\theta_{,zzz} - T^\theta_{,z} + T^\theta_{,\theta\theta z} + 2 u_{,\theta zz}\right)\nonumber\\ &\qquad\quad+\frac{\eta \epsilon}{r^3} \left(r^2 v^\theta_{,zz} + \partial_t u_{,\theta} + v^\theta_{,\theta\theta}\right)=0, \end{align}

(3.11)\begin{align} &\boldsymbol{e}_z:\quad\frac{\eta \epsilon}{r^2} \left(r^2 v^z_{,zz}-r \partial_t u_{,z} + v^z_{,\theta \theta}\right)+\zeta_{,z} = 0. \end{align}

The normal components of force due to the liquid crystal free energy are more complicated to calculate, however. Nevertheless, at $O(\epsilon )$, only the third term of (2.25) contributes, giving

(3.12)\begin{equation} \left(b^\alpha{}_\gamma T^\gamma T^\beta\right)_{;\beta\alpha} = \epsilon\left[u_{,zzzz}-\frac{1}{r^2}\,T^\theta_{,z\theta} +\frac{1}{r^2}\,u_{,zz\theta\theta}\right]+O(\epsilon^2), \end{equation}

which gives the shape equation

(3.13)$$\begin{gather} \frac{\mathcal{K} \epsilon}{r^2} \left(r^2 u_{,zzzz} - T^\theta_{z\theta} + u_{,zz\theta\theta}\right) -\frac{\zeta}{r^2} \left[ \epsilon\left(r^2 u_{,zz} + u_{,\theta\theta} + u\right)-r\right]\nonumber\\ {}+\frac{2 \eta \epsilon}{r^2} \left(\partial_t u + v^\theta_{,\theta}\right)=0. \end{gather}$$

Taking the case where $\epsilon =0$, we find a ground state with

(3.14)$$\begin{gather} \lambda=0, \end{gather}$$

(3.15)$$\begin{gather}\zeta = 0, \end{gather}$$

(3.16)$$\begin{gather}\boldsymbol{T}= \boldsymbol{e}_{z}, \end{gather}$$

as required. Note that as the bending energy term is anisotropic, it gives no characteristic length scale of the tube as compared to the case of a lipid bilayer, for example, where the ground-state radius of the tube is set by the ratio of isotropic bending energy to surface tension ($r_0=\sqrt {\kappa /2\zeta }$) (Zhong-Can & Helfrich Reference Zhong-Can and Helfrich1989; Nelson, Powers & Seifert Reference Nelson, Powers and Seifert1995; Gurin, Lebedev & Muratov Reference Gurin, Lebedev and Muratov1996; Fournier & Galatola Reference Fournier and Galatola2007; Powers Reference Powers2010; Boedec, Jaeger & Leonetti Reference Boedec, Jaeger and Leonetti2014).

We now take $\zeta =\epsilon \zeta$, $\lambda =\epsilon \lambda$, and take the Fourier transform in space $\bar {f}_{q,m}=(1/4{\rm \pi} ^2)\int \,\mathrm {d}z\,\mathrm {d}\theta \,f(z,\theta )\exp ({-\textrm {i} m \theta -\textrm {i} q z})$. Using each of (3.5), (3.6), (3.10), (3.9) and (3.13) leaves simply two dynamical equations for $\bar {u}_{qm}$ and $\bar {T}^\theta _{qm}$, which are given in dimensionless form as

(3.17) \begin{equation} \partial_{\tilde{t}} \left(\begin{matrix} \bar{T}^\theta_{qm}\\ \bar{u}_{qm} \end{matrix}\right)=\boldsymbol{\mathsf{A}} \cdot\left(\begin{matrix} \bar{T}^\theta_{qm}\\ \bar{u}_{qm} \end{matrix}\right) = \left(\begin{matrix} {\mathsf{A}}_{11} & {\mathsf{A}}_{12}\\ {\mathsf{A}}_{21} & {\mathsf{A}}_{22}\end{matrix}\right)\cdot\left(\begin{matrix} \bar{T}^\theta_{qm}\\ \bar{u}_{qm} \end{matrix}\right), \end{equation}

where componentwise,

(3.18) \begin{align} \left. \begin{aligned} {\mathsf{A}}_{11} = \frac{m^2 \left(-\left(Q^2 \left(4 \bar\eta +\nu ^2+1\right)+1\right)\right)-Q^2 \left(Q^2+1\right) \left(4\bar\eta +(\nu -1)^2\right)}{4 Q^2},\cr {\mathsf{A}}_{12} = \frac{m \left(m^4-2 m^2 \left((\nu -1) Q^2+1\right)+(1-2 \nu ) Q^4-2 Q^2 \left(4 \bar{\eta} +(\nu -1)^2\right)\right)}{4 Q},\cr {\mathsf{A}}_{21} = \frac{m \left({-}2 \nu Q^2 \left(m^2+Q^2+1\right)+m^2+Q^2\right)}{4 Q^3},\cr {\mathsf{A}}_{22}={-}\frac{m^6+m^4 \left(3 Q^2-2\right)+m^2 Q^2 \left(4 \nu +3 Q^2-2\right)+Q^6}{4 Q^2}, \end{aligned} \right\} \end{align}

where time has been normalised by $t=\tau \tilde {t}$ with $\tau =\eta r^2/\mathcal {K}$, $\bar {\eta }=\eta /\gamma$, $Q=rq$, and $\bar {u}$ has been non-dimensionalised by $r$. Note that the off-diagonal elements couple linearly in $mQ$. That is, the sign of $\bar {T}^{\theta }_{qm}$ selects a handedness of coupling to the shape deformations.

In contrast to the case of a flat membrane, this shows that in highly curved environments, the relaxation of shape and orientational ordering cannot be decoupled even in the case when the ground state has no intrinsic curvature. This is true for all the modes on the tube, with the exception of the $m=0$ case where the ordering and shape decouple. In figure 2, we plot the negatives of the eigenvalues $\{\varLambda _1,\varLambda _2\}$ of the matrix $\boldsymbol{\mathsf{A}}$ as functions of wavenumber $Q$ for the first $m=0,1,2,3$ modes. These can be interpreted as the passive relaxation modes that return perturbations back to the ground state.

Figure 2. Relaxation rates for a perturbed liquid-crystal tube. Here, $\varLambda _i$ are the two eigenvalues of the growth matrix given in (3.17), that is, the growth rate of the spatial Fourier modes with time non-dimensionalised by $\tau =\eta r^2/\mathcal {K}$, where $\eta$ is the 2-D viscosity, $r$ is the tube radius, and $\mathcal {K}$ is the liquid crystalline elastic modulus. These are plotted against dimensionless wavenumber $Q=q r$ ($z$ direction Fourier variable) for the first four integer modes in $m$ (the $\theta$ Fourier variable). The dimensionless viscosity is $\bar \eta =1$, and the spin connection is $\nu =-1.5$.

3.2. Active nematic liquid crystals

In this subsection, we will consider the effects of a non-zero active stress on the stability of nematic tubes, and on the morphology of a surface with a $+1$ defect located at the origin. In particular, we will examine the effects of changing the sign of the active stress (i.e. contractile versus extensile) on the morphology and surface instabilities.

3.2.1. Active ruffling and pearling instabilities in tubes

Using our equation for the active stress, we can find that the normal part of the active force for a tube is given by

(3.19)\begin{equation} \boldsymbol{f}^a\boldsymbol{\cdot}\hat{\boldsymbol{n}} = \sigma^{{(a)}\alpha\beta}b_{\alpha\beta}= \frac{\varSigma \left[\epsilon\left(r^2 u_{,zz}-u_{,\theta\theta} - u\right)+r\right]}{2 r^2}, \end{equation}

and the tangential parts by

(3.20)\begin{equation} (\mathbb{I}-\hat{\boldsymbol{n}}\otimes\hat{\boldsymbol{n}})\boldsymbol{\cdot}\boldsymbol{f}^a=\sigma^{{(a)}\alpha\beta}{}_{;\beta}\,\boldsymbol{e}_\alpha = \frac{\varSigma \epsilon T^\theta_{,z}}{r}\,\boldsymbol{e}_\theta + \frac{\varSigma \epsilon \left(2 r T^z_{,z} + T^\theta_{,\theta} +u_{,z}\right)}{r}\,\boldsymbol{e}_z. \end{equation}

where $\mathbb{I}$ is the identity matrix in $\mathbb{R}^3$.

Following a similar procedure as before, we now find a different ground-state surface tension at $\epsilon =0$ for the unperturbed tube of $\zeta = \varSigma /2 +O(\epsilon )$. Now following the same procedure as before, of Fourier transforming the equations around this ground state, we can find the dynamical equations for shape and alignment. Here, we will consider only the axisymmetric deformations ($m=0$) for simplicity. This leads to the following equations for $\bar {u}_{q}$ and $\bar {T}^\theta _{q}$:

(3.21)$$\begin{gather} \partial_{\tilde{t}}\bar{T}^{\theta}_{q}={-}\frac{(Q^2+1) \left((\nu -1)^2+4 \bar\eta \right)-2 (\nu -1) \bar\varSigma }{4}\,\bar{T}^\theta_q, \end{gather}$$

(3.22)$$\begin{gather}\partial_{\tilde{t}}\bar{u}_{q} ={-}\frac{Q^4+(Q^2-1) \bar{\varSigma} }{4}\,\bar{u}_q, \end{gather}$$

where $t=\tau \tilde {t}$, with $\tau =\eta r^2/\mathcal {K}$, $\bar {\eta }=\eta /\gamma$, $\bar {\varSigma }=\varSigma r^2/\mathcal {K}$, $Q=rq$, and $\bar {u}$ has been non-dimensionalised by $r$. Note that $\bar \varSigma$ is essentially a ratio of surface forces to bending energy and can thus be viewed as an ‘active Föppl–von Kármán number’ that we discussed earlier.

First, we will consider the second equation describing the shape dynamics. The growth rate for the shape perturbation, $\bar {u}_q$, is given by

(3.23)\begin{equation} \varLambda_u(Q,\bar{\varSigma})={-}\frac{Q^4+\left(Q^2-1\right) \bar{\varSigma}}{4}. \end{equation}

This can be driven unstable for

(3.24)\begin{equation} \bar\varSigma\leq \frac{-Q^4}{Q^2-1}. \end{equation}

For the case of an extensile active stress ($\varSigma <0$), this leads to an instability of a finite range of $Q>1$ peaked around a specific wavelength; see figure 3. The fastest growing mode of this instability is found at $Q^\star =\sqrt {-\bar \varSigma /2}$, and as this occurs at $Q>1$, leads to ruffling-like undulations in the membrane. This makes sense heuristically as the extensile stresses are pointing along the tube thus generating forces that act to increase the surface area of the tube, hence the ruffling.

Figure 3. Extensile instability causes a ruffling-like deformation in the shape of the tube. Here, $\varLambda _u(Q,\bar {\varSigma })$ is the growth rate of the spatial Fourier modes, with time non-dimensionalised by $\tau =\eta r^2/\mathcal {K}$, where $\eta$ is the 2-D viscosity, $r$ is the tube radius, and $\mathcal {K}$ is the liquid crystalline elastic modulus. The active stress is non-dimensionalised as $\bar \varSigma =\varSigma r^2/\mathcal {K}$. Inset shows a representative example of such an instability.

In the case of contractile stresses ($\varSigma >0$), we find an instability with fastest growing mode in the $Q\to 0$ limit; see figure 4(a). This instability is analogous to a classical Rayleigh–Plateau instability (Rayleigh Reference Rayleigh1892; Tomotika Reference Tomotika1935) and similar to the pearling instabilities seen in fluid membrane tubes (Gurin et al. Reference Gurin, Lebedev and Muratov1996; Boedec et al. Reference Boedec, Jaeger and Leonetti2014; Narsimhan, Spann & Shaqfeh Reference Narsimhan, Spann and Shaqfeh2015). In reality, such an instability does not occur at infinite length scales as the dissipative time scale in the ambient media damps the longer-wavelength modes. We can approximate the bulk dissipation by modifying our shape dynamics equation as (Nelson et al. Reference Nelson, Powers and Seifert1995; Gurin et al. Reference Gurin, Lebedev and Muratov1996; Powers Reference Powers2010; Al-Izzi et al. Reference Al-Izzi, Sens, Turner and Komura2020b; Boedec et al. Reference Boedec, Jaeger and Leonetti2014)

(3.25)\begin{equation} \varLambda_u(Q,\bar{\varSigma})={-}\frac{Q^4+\left(Q^2-1\right) \bar{\varSigma}}{4+8 Q^{{-}2}\tilde{L}_{{SD}}^{{-}1}}, \end{equation}

where $\tilde {L}_{{SD}} = \eta /(\eta _b r)$ is the dimensionless Saffman–Delbrück length, and $\eta _b$ is the ambient viscosity (Saffman Reference Saffman1975; Saffman & Delbrück Reference Saffman and Delbrück1975). Here, this leads to a long-wavelength instability as the contractile stress tries to minimise the surface to volume ratio of the tube. This gives a finite wavelength to the fastest growing mode set by the Saffman–Delbrück length; see figure 4(b). This instability is similar to the isotropic contractile instabilities found previously in fluid membranes (Mietke, Julicher & Sbalzarini Reference Mietke, Jülicher and Sbalzarini2019).

Figure 4. Contractile instability causes a pearling-like deformation in the shape of the tube. Here, $\varLambda _u(Q,\bar {\varSigma })$ is the growth rate of the spatial Fourier modes, with time non-dimensionalised by $\tau =\eta r^2/\mathcal {K}$, where $\eta$ is the 2-D viscosity, $r$ is the tube radius, and $\mathcal {K}$ is the liquid crystalline elastic modulus. The active stress is non-dimensionalised as $\bar \varSigma =\varSigma r^2/\mathcal {K}$. (a) Instability in the absence of an ambient fluid. (b) Instability when an approximation to ambient hydrodynamics in the long-wavelength limit is used. Here, $\tilde {L}_{{SD}}= {\eta }/{\eta _b r}$ is the dimensionless Saffman–Delbrück length given by the ratio of 2-D to three-dimensional viscosity ($\eta _b$) non-dimensionalised by tube radius. Inset shows an example growth of a surface undulation for the case with bulk fluid.

We now turn our attention to the growth rate of the director deviations, which is given by

(3.26)\begin{equation} \varLambda_T(Q,\bar{\varSigma},\bar{\eta},\nu)={-}\frac{\left(Q^2+1\right) \left(4 \bar\eta +(\nu -1)^2\right)-2 (\nu -1) \bar\varSigma}{4}. \end{equation}

This can also have a change in sign when in the extensile or contractile regime, depending on the sign of $\nu$ ($\nu <-1$ corresponds to rod-like nematics, and $\nu >1$ to disk-like). The stability condition is given by

(3.27)\begin{equation} \bar\varSigma\geq \frac{(Q^2+1)\left(4\bar\eta+(\nu-1)^2\right)}{2(\nu-1)}, \end{equation}

which leads to a spontaneous bend instability in the longest wavelengths of the tube; see figure 5. This instability is well known in flat geometry in the case of extensile stress and rod-like nematics ($\nu <-1$), where it is known to lead to active turbulence (Marchetti et al. Reference Marchetti, Joanny, Ramaswamy, Liverpool, Prost, Rao and Simha2013).

Figure 5. (a) Extensile bend instability in texture on the surface of the tube with rod-like nematics. (b) Contractile bend instability in texture on the surface of the tube with disk-like nematics. Here, $\varLambda _T(Q,\bar {\varSigma })$ is the growth rate of the spacial Fourier modes, with time non-dimensionalised by $\tau =\eta r^2/\mathcal {K}$, where $\eta$ is the 2-D viscosity, $r$ is the tube radius, and $\mathcal {K}$ is the liquid crystalline elastic modulus. The active stress is non-dimensionalised as $\bar \varSigma =\varSigma r^2/\mathcal {K}$. The dimensionless viscosity is $\bar \eta =1$, and the spin connection is (a) $\nu =-1.5$, and (b) $\nu =1.5$. Insets show an example of the bend instability at the longest wavelength of the tube.

In the limit $Q\to 0$, this criterion becomes $\bar \varSigma \geq (4\bar \eta +(\nu -1)^2)/2(\nu -1)$, in comparison to the flat case where the instability is thresholdless (Marchetti et al. Reference Marchetti, Joanny, Ramaswamy, Liverpool, Prost, Rao and Simha2013). This is due directly to the correct use of the surface derivative as the texture needs sufficient active stress to overcome the elastic stress the bend induces by forcing the texture to bend around the tube. Because of this threshold, it is interesting to note that the shape instability precedes the bend instability, suggesting that it should be possible to generate shape changes before one sees active turbulence. The interplay between this instability and the shape equations beyond linear order is likely a very rich topic as this director instability would likely couple to helical deformation modes leading to a spontaneous chiral symmetry breaking. However, such a topic is beyond the scope of the analysis in this paper.

3.2.2. Active topological defect deformations

Here, we consider a planar membrane in polar Monge coordinates, where the surface is parametrised by the Cartesian vector $\boldsymbol {R}=(r\cos \theta,r\sin \theta,\epsilon z)$. For simplicity, we will consider only axisymmetric steady states, thus $z(r,\theta )=z(r)$. The metric and second fundamental form are given as

(3.28)\begin{gather} {}[g_{\alpha\beta}] = \left( \begin{array}{@{}cc@{}} 1 & 0 \\ 0 & r^2 \end{array} \right) +O(\epsilon^2), \end{gather}

(3.29)\begin{gather} {}[b_{\alpha\beta}] = \left( \begin{array}{@{}cc@{}} \epsilon z_{,rr} & 0 \\ 0 & r \epsilon z_{,\theta} \end{array} \right) + O(\epsilon^2), \end{gather}

up to first order in $\epsilon$. Thus we find $H= ({\epsilon }/{2})(({1}/{r})z_{,r}+z_{,rr})+O(\epsilon ^2)$ and $K=0+O(\epsilon ^2)$. We take the polarisation texture to be $\boldsymbol {T}=\cos \phi \,\boldsymbol {e}_{r} + ({1}/{r})\sin \phi \,\boldsymbol {e}_\theta$, where $\phi$ is a constant such that there is a $+1$ topological defect located at $r=0$. The continuity equation then yields $v^r=0$, so we are left to solve for $z(r,t)$, $\zeta (r,t)$, $\lambda (r,t)$ and $v^\theta (r,t)$.

At lowest order in $\epsilon$, the equations are independent of $\dot {z}$, so we look for steady states in the nematic order parameter. Our choice of $\boldsymbol {T}$ trivially satisfies the unit magnitude constraint. We now make use of the fact that at lowest order, the tangential force balance and polarisation equations are just given by their flat counterparts (that is, they have no $O(\epsilon )$ term), and the shape equation is given entirely by $O(\epsilon )$ terms (Kruse et al. Reference Kruse, Joanny, Jülicher, Prost and Sekimoto2005; Hoffmann et al. Reference Hoffmann, Carenza, Eckert and Giomi2022).

The polarisation equations become

(3.30)$$\begin{gather} (\nu \cos (2 \phi )-1) \left( v^\theta_{,r}-\frac{v^\theta}{r}\right)=0, \end{gather}$$

(3.31)$$\begin{gather}\nu \sin (\phi ) \cos (\phi ) \left( v^\theta_{,r}-\frac{v^\theta}{r}\right)+\lambda +\frac{\mathcal{K} }{\gamma r^2}=0, \end{gather}$$

which, assuming $\nu >1$ (disk-like nematics), gives

(3.32)$$\begin{gather} \phi ={\pm} \frac{1}{2}\arccos\left[\frac{1}{\nu}\right], \end{gather}$$

(3.33)$$\begin{gather}\lambda ={-} \frac{\mathcal{K}}{r^2 \gamma} - \nu\sin (\phi)\cos(\phi)\left[v^\theta_{,r}-\frac{v^\theta}{r}\right]. \end{gather}$$

Substituting the solution for $\phi$ into the tangential force balance gives

(3.34)$$\begin{gather} \zeta_{,r} + \frac{\mathcal{K} \nu }{r^3}+\frac{\varSigma }{\nu r} = 0, \end{gather}$$

(3.35)$$\begin{gather}r \left(\eta v^\theta_{,r} + \eta r v^\theta_{,rr} \pm \sqrt{1-\frac{1}{\nu ^2}}\,\varSigma \right)-\eta v^\theta = 0, \end{gather}$$

for which solutions are

(3.36)$$\begin{gather} \zeta = \frac{\mathcal{K}\nu}{2 r^2} - \frac{\varSigma}{\nu} \log\left[\frac{r}{r_{{max}}}\right]+\zeta_0, \end{gather}$$

(3.37)$$\begin{gather}v^\theta ={\mp} \frac{\varSigma}{2\eta} \sqrt{1-\frac{1}{\nu^2}}\,r \log\left[\frac{r}{r_{{max}}}\right], \end{gather}$$

where $\zeta _0 + \mathcal {K}\nu /2 r^2_{{max}}$ is the surface tension at $r=r_{{max}}$ (figure 6). This solution is identical to that found in Kruse et al. (Reference Kruse, Joanny, Jülicher, Prost and Sekimoto2004), and is one of the classical results of polar active gel theory in a flat geometry. Notice that the surface tension diverges as $r\to 0$ due to the presence of the defect. This divergence in shape can be avoided by introducing a short wavelength cut-off $r_c$, or by considering a more general Landau–de Gennes free energy that allows for a nematic isotropic phase transition (DeGennes & Prost Reference DeGennes and Prost1993).

Figure 6. (a) Defect surface tension $\zeta (r)$ as a function of radius $r$ about a $+1$ defect centred on $r=0$. Inset is defect texture. (b) Velocity field around the defect. Here, we have introduced a small cut-off radius at $r=0.03$, where lengths are non-dimensionalised by the radius of the boundary $r_{{max}}=1$. Surface tension is given by $\bar \zeta _0=\zeta _0 r_{{max}}^2/\mathcal {K}=1.0$, spin connection is $\nu =1.5$, $r_c=0.03$ and $c=-0.05$, where $c$ essentially defines the contact angle at the boundary. The active stress is non-dimensionalised by $\bar \varSigma = \sigma r_{{max}}/\mathcal {K}$. We plot for both positive (contractile) and negative (extensile) values of active stress to show promotion/suppression of the deformation in each case.

Substituting our solution for the surface tension into the normal force balance equation yields the dimensionless shape equation

(3.38)$$\begin{gather} \frac{z_{,r}}{r^3} \left(7-\nu \left(2 \bar\zeta_0 r^2+5\right)+2 r^2 \bar\varSigma \log (r)+r \left(\sqrt{\nu ^2-1}+r \bar\varSigma \right)\right)\nonumber\\ {}-\frac{z_{,rr}}{r^2} \left(r^2 (2 \bar\zeta_0 \nu + \bar\varSigma )-2 r^2 \bar\varSigma \log (r)+4\right)\nonumber\\ +\,(\nu +1) \left( z_{,rrrr} + \frac{2}{r}\,z_{,rrr}\right)=0, \end{gather}$$

where lengths have been non-dimensionalised by $r_{{max}}$, and energies by $\mathcal {K}$. We can solve this equation numerically from the boundary ($r=1$ in our dimensionless units) up to a small cut-off $r_c$. We solve with zero mean curvature and derivative in curvature at $r=r_{{max}}=1$ – at lowest order, this is $z'(r_{{max}})=-z''(r_{{max}})=z'''(r_{{max}})=c=\textrm {Const.}$ – and fix $z(r_{{max}})=0$ at the boundary. The constant $c$ essentially prescribes the contact angle at the boundary. To solve this numerically, we make use of Mathematica (WolframResearch, Champaign, IL).

Solutions to this for varying values of active stress are plotted in figure 7 with $\bar \zeta _0=1.0$, $\nu =1.5$, $r_c=0.03$ and $c=-0.05$. For $\bar \varSigma =0$, we find a convex cusp solution where the stress from the defect causes the surface to buckle. Note that this is qualitatively similar to the full nonlinear pseudo-sphere cusp found in (Frank & Kardar Reference Frank and Kardar2008), although the quantitative details of our shape differ significantly due to our inclusion of the extrinsic curvature term in the Frank free energy, and thus the anisotropic bending energies in our shape equation (fourth- and third-order terms) making a full analytical nonlinear solution of our equations impossible. For the extensile case, the protrusion decreases in height. The reason for this is due directly to our inclusion of the extrinsic terms in the Frank free energy and the convex shape of the passive solution. Extensile stresses lead to active forces that act to increase bend, in this case bend due to the deformation along the surface, thus pushing the surface down and suppressing the deformation. The opposite is true in the contractile case, where the active forces act to flatten bend, thus pushing the surface upwards. A simple example of these types of effect was discussed in the context of an active fingering instability in Alert (Reference Alert2022). In addition, we show that as the spin connection $\nu$ (or liquid crystalline Scriven–Love number) is increased, the spiral defect texture goes from aster-like to vortex-like as $\nu$ goes to $\infty$. As this happens, the bend along the surface goes to zero, and the height of the defect, $z(r_{c})$, converges in the extensile, contractile and passive cases, as expected.

Figure 7. (a) Surface heights $z(r)$ as a function of radius $r$ about a $+1$ defect centred on $r=0$. Here, we have introduced a small cut-off radius at $r=0.03$, where lengths are non-dimensionalised by the radius of the boundary $r_{{max}}=1$. Surface tension is given by $\bar \zeta _0=\zeta _0 r_{{max}}^2/\mathcal {K}=1.0$, spin connection is $\nu =1.5$, $r_c=0.03$ and $c=-0.05$, where $c$ essentially defines the contact angle at the boundary. The active stress is non-dimensionalised by $\bar \varSigma = \sigma r_{{max}}/\mathcal {K}$. We plot for both positive (contractile) and negative (extensile) values of active stress to show promotion/suppression of the deformation in each case. Inset is an example surface with the texture shown on top. (b) Maximum defect height $z(r_c=0.03)$ plotted against spin connection $\nu$ (or liquid crystalline Scriven–Love number) in the passive, contractile and extensile cases.

Note that our results here are different from those seen in Hoffmann et al. (Reference Hoffmann, Carenza, Eckert and Giomi2022), where the effects of contractile and extensile forces are switched. This is straightforward to understand since Hoffmann et al. (Reference Hoffmann, Carenza, Eckert and Giomi2022) did not consider the extrinsic coupling of the director, and in addition included an isotropic bending energy. This causes the passive shape to be regularised around the defect and to be concave near the defect where the active forces are strongest, leading to a switch in sign of the normal force due to activity when compared to our convex passive shape. Such subtle interplay between passive and active stresses suggests a rich phenomenology of possible nonlinear solutions to such systems, and that it may be possible for both extensile and contractile systems to achieve a similar morphology by the varying of passive parameters.

It is interesting to note that the contractile protrusions we find are similar to those seen during Hydra morphogenesis in the underlying contractile actomyosin network (Maroudas-Sacks et al. Reference Maroudas-Sacks, Garion, Shani-Zerbib, Livshits, Braun and Kinneret2021; Vafa & Mahadevan Reference Vafa and Mahadevan2022). There, the $+1$ defects form the basis of protrusions that develop into ‘limbs’.