Self-organized shape dynamics of active surfaces.

Mechanochemical processes in thin biological structures, such as the cellular cortex or epithelial sheets, play a key role during the morphogenesis of cells and tissues. In particular, they are responsible for the dynamical organization of active stresses that lead to flows and deformations of the material. Consequently, advective transport redistributes force-generating molecules and thereby contributes to a complex mechanochemical feedback loop. It has been shown in fixed geometries that this mechanism enables patterning, but the interplay of these processes with shape changes of the material remains to be explored. In this work, we study the fully self-organized shape dynamics using the theory of active fluids on deforming surfaces and develop a numerical approach to solve the corresponding force and torque balance equations. We describe the spontaneous generation of nontrivial surface shapes, shape oscillations, and directed surface flows that resemble peristaltic waves from self-organized, mechanochemical processes on the deforming surface. Our approach provides opportunities to explore the dynamics of self-organized active surfaces and can help to understand the role of shape as an integral element of the mechanochemical organization of morphogenetic processes.

Morphogenesis is the generation of patterns and shapes by dynamic processes. In the biological context, cells and tissues are shaped during developmental processes that give rise to organisms starting from a fertilized egg. Biological morphogenesis is fundamentally mechanochemical; i.e., it relies on an interplay between chemical signals and active mechanics (1–3). Active mechanical processes in cells are generated in the cytoskeleton, a gel-like network of protein filaments. Motor proteins, such as myosin, interact with cytoskeletal filaments and generate movement and forces driven by the hydrolysis of ATP (4), thus making the cytoskeleton an active material. The capacity of living matter to generate spontaneous motion, flows, and material deformations as a result of molecular force-generating processes can be captured by the theory of active matter (5, 6).

The cellular actomyosin cortex is a thin layer of an active material near the cell membrane, which governs cell mechanics and cell shape changes. For example, during cell division a contractile ring forms in the cell cortex that drives cortical flows and constricts the cell to create two daughter cells (7–10). Active cortical contractions also play a role in tissue shape changes (11, 12). For example, during gastrulation it has been suggested that contraction of the apical cell cortex drives the invagination of a 2D cell layer, called the epithelium (13). These examples show that active mechanical processes on surfaces play a key role in driving morphogenesis.

The emergent character of morphogenetic processes naturally arises from a feedback loop in which chemically regulated active stresses induce material flows and deformations that in turn affect the chemical regulators (14). This raises the fundamental question of how active stresses are dynamically organized during shape changes of cells and tissues. In the absence of shape changes, it has already been demonstrated that advection of regulators of active stress gives rise to mechanochemical self-organization and patterning (14–16). The general process of morphogenesis introduces additional feedback via changes of the underlying geometry that enables the generation of shapes by self-organizing processes. Here, we present a general framework to describe the dynamics of shapes that arise from the self-organization of mechanochemical processes on an active surface and develop a numerical technique to capture the corresponding shape changes of surfaces with axial symmetry.

Our formulation can be seen as a generalization of previous work, in which shapes are determined from a minimization of energy functionals (17–19) or from a deformations dynamic that is guided by phenomenological rules (20–22). Here, we derive the deformation dynamics and the steady states from the general force and torque balance on curved surfaces. This allows including constitutive relations that are associated with energy functionals of the surface shape (23), as well as constitutive relations that describe active elastic or viscous materials (24, 25). In particular, we consider in this work an active thin film description, which has been extensively applied as a generic model to account for the dynamic behavior of the cellular cortex and tissue layers (7, 16, 26, 27). In our description, active stresses in the material can depend on the concentration of a stress-regulating molecular species that is dynamically changing in response to flows on and deformations of the surface. We study this system on surfaces with spherical and tubular geometries, which represent a broad class of biological structures including the cellular cortex (7, 9, 17, 24, 26, 27) and various types of biological tissue (12, 16, 19, 28).

To solve the resulting equations and the shape dynamics numerically, we use an integral representation of the deforming surface and develop a dynamic coordinate transformation. We show that a mechanochemical instability of a spherical active surface leads to spontaneous polarization and the formation of deformed stationary shapes. We point out that surface deformations, even in the absence of this instability, lead to an accumulation of stress regulators in regions of high curvature due to geometric effects. Furthermore, we show that the self-organization of active stress regulators can spontaneously constrict tubular surfaces and can lead to mechanochemical shape oscillations. Finally, we describe the emergent formation of directed flows on tubular surfaces with bending rigidity that resemble a peristaltic wave resulting from the coupling between self-organized active stresses and surface geometry.

Geometry and Mechanics of Self-Organized Active Surfaces

We consider a time-dependent surface that is represented by a parameterization . The parameters and denote the surface coordinates and time, respectively. Basis vectors in the tangent space of are given by (Fig. 1), where denote partial derivatives with respect to the surface coordinates. The metric tensor is defined by and , where we use the Einstein summation convention. With the surface normal , vector fields on the surface can be decomposed into tangential and normal parts , where and . The surface shape is characterized by the curvature tensor . The mean curvature of the surface is defined by .

Fig. 1.

Description of curved surfaces. (A) Parameterization of a curved surface embedded in . The tension tensor and moment tensor are used to describe forces and torques in the surface acting on a line element . (B) Representation of an axisymmetric, deforming surface with arc-length parameterization map . Shown is a schematic set of collocation points used for the numerical discretization, mapped onto the inhomogeneous grid in physical space via the coordinate transformation . The meridional curvature encodes the shape of axisymmetric surfaces. The vectors denote the normalized standard basis in cylindrical coordinates, is the tangent angle, and .

Force and Torque Balance Equations of Curved Surfaces.

To describe deformations of as a result of mechanical stresses within and onto the surface, surface configurations should obey force and torque balance equations that we briefly introduce in the following. Consider a line element of length across the surface with in-plane normal (Fig. 1). The forces acting on this line element can be written as (23). This defines the in-plane tension tensor and the out-of-plane shear tension by . Torques acting on the line element can be treated in the same way. Additionally, external forces may act on the surface. With these definitions, the force and torque balance equations read (23–25)where denotes the covariant derivative and denotes the covariant Levi–Civita tensor (). Eqs. and describe the force balance in the directions tangential and normal to the surface, respectively, and we have omitted inertial effects. Eq. determines how bending moments in the surface are balanced by out-of-plane shear tension , and Eq. determines their coupling to the antisymmetric part of the tension tensor . Constitutive relations of the material provide expressions for the tension and moment tensors (25), which are introduced in .

Constitutive Relations of Active Fluid Surfaces.

To define the active thin film description studied in this work, we first introduce the equilibrium properties of the surface. We consider here a surface with constant surface tension and bending rigidity , described by the Helfrich free energy (29)Bending rigidity has been implicated as a key element in determining cell and tissue shapes (3, 19, 21). The tension and moment tensors of a surface described by Eq. read (23)Spontaneous curvature and Gaussian bending rigidity lead to additional terms in Eqs. –, which we do not consider here for simplicity (25). Note that the normal force balance Eq. yields for and the well-known shape equation of membranes (23).

To characterize flows and deformations of the thin film, we use the symmetric part of the in-plane strain rate tensor , which reads (25)The component denotes the in-plane flows, describes deformations of the surface, and is the flow velocity of surface elements. The deformation velocity contributes to strains that are generated when deforming curved surface regions, such that curvature becomes a natural element of the constitutive equations of a deforming surface. The deviatoric tension of the isotropic active thin film is given by (25, 27)Here, and are the shear and bulk viscosities. The viscous limit is justified whenever the turnover time scales of the substructures that constitute the surface are smaller than the time scales of deformations and in-plane flows (8, 30, 31). An effective in-plane compressibility can for example arise from surface height changes (8, 32) or from an exchange of surface material with the environment, which we do not explicitly take into account. The last term with describes an isotropic, contractile active tension in the surface. The constitutive relations can be obtained usingwith the contributions to the tension and moment tensors given by Eqs. , , and . Additionally, we consider and note that can be computed using the moment balance Eq. . Furthermore, we assume that the volume enclosed by the surface is conserved, which defines a pressure that enters the normal force balance Eq. as ().

Chemical Regulation of Active Tension.

We consider an active tension amplitude of the form , where is the contractility and is a concentration field that regulates the strength of the local active tension. Following previous works (7, 14, 15), this regulation is represented by a saturating Hill functionThe concentration field is changing over time due to advective transport, reactions, and diffusion and therefore plays a crucial role for the self-organization of morphogenetic processes. These processes are captured by the dynamic equation for the concentration field (25, 33)The first term on the right-hand side describes advection due to in-plane flows. The second term describes dilution and accumulation due to local surface expansion () or compression (), respectively, that occurs during deformations of surface regions with nonvanishing mean curvature . Isotropic in-plane diffusion with diffusion constant is incorporated using the Laplace–Beltrami operator . The last term describes the turnover with a rate and preferred concentration of the stress regulator .

Dynamic Representation of Deforming Axisymmetric Surfaces

For a given deformation velocity , the dynamics of the surface represented by the Lagrangian parameterization obeyswhere denotes the total time derivative. In the following, we consider deformations of axisymmetric surfaces and develop a framework to numerically solve Eq. for a given deformation velocity .

An arbitrary axisymmetric surface can be represented explicitly using an arc-length parameterizationwhere is the normalized standard basis with in cylindrical coordinates, is the azimuthal angle, and is the arc-length parameter of the meridional outline (Fig. 1). In the following, we use and explicitly as tensor indexes. The nonvanishing components of the curvature tensor are the azimuthal curvature and the meridional curvaturewhere is the tangent angle defined by . This definition implies that the full shape information about is encoded in the meridional curvature , together with a point . Indeed, for given , we can compute the tangent angle from Eq. asand reconstruct , parameterized as in Eq. , viaConsequently, the shape of the surface during deformations, and therefore the solution to Eq. , can be found by determining the time evolution of the meridional curvature and of the values of , , and at .

An arc-length parameterization simplifies the parametric form of covariant equations and the shape reconstruction via Eqs. –. However, the time dependence of the domain makes it difficult to evaluate the total time derivative in the equation of the shape dynamics, Eq. , and renders the arc-length parameterization impractical for the numerical treatment of differential equations on a deforming surface. We therefore introduce additionally a Eulerian parameterization of given bywhere is a parameter on an interval that remains fixed during surface deformations. For convenience, we chose here , such that . The azimuthal angle is defined as in Eq. . Note that and represent the same surface , but for the Eulerian parameterization, the shape dynamics Eq. can be written asThe surface coordinates and are then related by a time-dependent coordinate transformation that is defined byTo obtain a dynamic equation for , we first note that we have and thus . Furthermore, it follows from Eq. and the definition of the metric tensor that for an Eulerian parameterization () , which impliesHere, we have used . For our choice of the fixed interval where is defined, is uniquely determined as the solution of Eq. with initial condition .

To reconstruct the deforming surface via Eqs. – as a function of time, we need the meridional curvature that obeys ()The first term captures curvature changes due to stretching of the meridional outline, and the second term describes the outline’s bending that occurs in the case of an inhomogeneous deformation velocity . Altogether, Eqs. – and – provide a framework to solve Eq. for a given deformation velocity on axisymmetric surfaces.

Self-Organization of Active Surface Deformations

We now combine the dynamics of surfaces with the physics of an active fluid film. For a given surface shape and distribution of active tension, the instantaneous deformation velocity and the meridional in-plane flow can be obtained by solving hydrodynamic equations of the active fluid surface. The hydrodynamic equations result from combining the force and torque balance Eqs. – with the constitutive relations Eqs. and (). Note that these equations determine the velocities up to a constant velocity vector because the system is Galilei invariant. The dynamics of the shape together with dynamics of the stress regulator are then obtained by solving Eqs. and , using at all times the instantaneously determined velocity fields. This combination of geometry and active hydrodynamics captures the mechanochemical feedback mediated by the chemical regulator, where the latter organizes active stress patterns and depends itself on material flows and on changes of the geometry.

We now study self-organized surface deformations of axisymmetric spherical and tubular surfaces. To characterize the strength of the mechanochemical feedback, we introduce a dimensionless contractility parameterIn the following, we first perform a linear stability analysis of homogeneous states to determine parameter regimes where the surface is unstable. We then study the full nonlinear surface dynamics using the framework introduced in the preceding sections and detailed further in .

Self-Organized Shape Dynamics of Spherical Surfaces.

We first consider axisymmetric surfaces of spherical topology and analyze the linear stability of the homogeneous steady state with and on a sphere with radius . We expand small perturbations of this state as , , and . Here, , , and denote axisymmetric scalar and vector spherical harmonics (34). Expanding the force balance Eqs. – and the dynamic Eq. for the concentration field to linear order (), we find that for increasing contractility parameter the mode has the largest growth rate and is the first mode that becomes unstable at , whereThis instability is independent of the bending rigidity .

The growing mode characterizes a polar asymmetry. Using a randomly perturbed concentration field as initial condition in our numerical approach, we find that for the instability leads to the spontaneous formation of a single patch of the stress regulator and an asymmetric surface shape (Fig. 2 and Movie S1). In the final steady state the advective influx of stress regulator into the contractile region is balanced by a diffusive outflux away from it. The resulting inhomogeneous tension across the surface leads to an oblate shape with broken mirror symmetry with respect to the axis and thus spatially varying curvature (Fig. 2, Right).

Fig. 2.

Shape dynamics of spherical surfaces with concentration-dependent active tension. (A) Schematic stability diagram of the system. The critical contractility is independent of the bending rigidity . (B) Mechanochemically unstable surface with and . (Left) After a small concentration perturbation on a sphere (, not to scale), a deformed steady-state shape emerges with a localized patch of stress regulator (Movie S1). Red arrows denote the in-plane flow field . (Right) Profiles of the concentration , in-plane flow (, ), and principle curvatures , of the steady-state surface. In the steady state, the diffusive outflux away from the contractile patch is balanced by an advective influx. The resulting tension across the surface is inhomogeneous, which leads to a deformed steady-state shape. (C) Relaxation dynamics of a mechanochemically stable surface with and , starting with a spheroidal shape with eccentricity 0.75 (Movie S2). Inhomogeneities in the initial mean curvature () lead to transient inhomogeneities in the concentration field () due to a deformation-induced local expansion and compression of the surface, before the stable steady state of a sphere is reached (, not to scale). Surface flows in B and C are shown in the reference frame where the pole velocities satisfy (). The parameters used in these simulations are given in .

For the homogeneous concentration on the sphere represents a stable solution. This is revealed by the relaxation of a deformed sphere toward this stable steady state. Such a relaxation is depicted in Fig. 2, where we use a spheroidal surface with a homogeneous concentration as initial condition. During the relaxation process, we observe the transient formation of concentration maxima at the poles (Movie S2). These maxima appear as a consequence of the large mean curvature at those locations, which leads to a locally increased surface compression during deformations (Eq. ).

Self-Organized Shape Dynamics of Tubular Surfaces.

We now study the self-organization of an active surface with a tubular geometry. We analyze the linear stability of the homogeneous steady state with and on a cylinder with radius and length and consider periodic boundary conditions. We expand small axisymmetric perturbations of this state as , , and with . From the corresponding expansion of the force balance Eqs. – and the dynamic Eq. for the concentration field to linear order, we derive the stability diagram as a function of the contractility parameter and the aspect ratio (Fig. 3 and ). For the cylinder surface is unstable in the blue-shaded region. For withthe instability occurs when the aspect ratio increases beyond a critical value that equals and does not depend on . This instability corresponds to the Plateau–Rayleigh instability and is not driven by the self-organization of the stress regulator. The case and exactly represents the classical Plateau–Rayleigh scenario (35). For instabilities are driven by the mechanochemical self-organization of the fluid surface. In this regime the aspect ratio at which the cylinder surface becomes unstable is smaller than and decreases for increasing contractility parameter . The linear stability analysis reveals furthermore a region where eigenvalues are complex (dark blue-shaded area in Fig. 3), which indicates oscillatory behavior at the instability. These characteristics of the stability diagram remain qualitatively unchanged when the bending rigidity is finite. The instability line for is shown as a red dashed line in Fig. 3. The red-shaded region indicates complex eigenvalues. Note that the value of in general depends on ().

Fig. 3.

Shape dynamics of tubular active fluid surfaces with concentration-dependent active tension. (A) Representative stability diagram of the system. Stable and unstable regions are separated by the blue curve for and the red dashed line for , where . For the dark blue-shaded region indicates parameter regimes where eigenvalues of the Jacobian are complex. The critical contractility (Eq. ) is described in the main text. (B) Concentration and in-plane flow (red arrows) during the spontaneous formation of a contractile ring (Movie S3) for parameters indicated with a blue circle in the stability diagram (, not to scale). This mechanism can constrict surfaces with aspect ratios below the Rayleigh threshold of . (C) Concentration and in-plane flow over one oscillation period for parameters indicated with a blue cross in the stability diagram (Movie S4). Oscillations result from the interplay between geometric stability of cylinder surfaces with and the mechanochemical instability of the active fluid film. Surface flows in B and C are shown in the reference frame where (). (D) Steady state of directed surface flows relative to a constricted shape ( not to scale). Shown are concentration (color code as in B), in-plane flow (), and principle curvatures , . The surface flow is shown in the reference frame where the constriction does not move (). In the reference frame where , this steady state resembles peristaltic motion (Movie S5). The parameters used in these simulations are given in .

We also study the surface dynamics beyond the linear regime, using our numerical approach. For , , and aspect ratios inside the unstable region of the stability diagram (blue circle in Fig. 3), the cylinder surface constricts and generates a thin cylindrical neck region with decreasing radius (Fig. 3 and Movie S3). The numerical analysis indicates that this radius vanishes at finite time. The concentration of the stress regulator increases along the tubular neck. For parameters that correspond to complex eigenvalues in the linear stability diagram (blue cross in Fig. 3), the cylinder constricts and expands periodically (Fig. 3 and Movie S4) with increasing amplitude until the neck radius vanishes. For and , the instability also leads to a tubular neck with vanishing radius, however, not driven by contractility but by surface tension due to the Plateau–Rayleigh character of the instability in this parameter regime. In this case, the concentration of stress regulator increases in the neck predominantly because of a local reduction of surface area.

For , an unstable cylinder surface constricts and reaches at long times a finite neck radius for all values of the contractility parameter . For a high concentration of stress regulator additionally builds up in the neck region. The concentration pattern and surface flows that emerge at long times spontaneously break the mirror symmetry with respect to the axis (Fig. 3 and Movie S5). As a consequence, average flows directed along the axis occur in a reference frame, where the constriction does not move. This corresponds to a peristaltic contraction wave that propagates in a reference frame where the average flow vanishes. Our numerical results reveal that such propagating solutions emerge in all parameter regimes for which the cylinder surface is linearly unstable and . For no propagating solutions occur. To test whether the stabilization of a finite neck radius is sufficient for directed flows to occur via spontaneous symmetry breaking, we chose and consider an external force that also stabilizes small neck radii. In this setting, we again find the formation of directed flows relative to the constricted shape. This shows that the spontaneous symmetry breaking for is an intrinsic property of the self-organized active fluid film and its tight coupling to the surface geometry.

Discussion

We have introduced a simple but general model for the mechanochemical self-organization of surface geometry. Active stresses in the surface are regulated by a diffusible and advected molecular species. Gradients of active stress induce surface flows and shape changes, which in turn influence the distribution of the stress regulator. As a consequence, shape changes, shape oscillations, and spontaneous surface flows can be generated via mechanochemical instabilities. In contrast to mechanochemical instabilities that have been discussed in fixed geometries (14–16), the phenomena described in the present paper give rise to shape changes and depend themselves crucially on the shape changes that occur.

To solve the dynamic equations for the shape, the flows, and the concentration fields, we have developed a numerical approach based on an integral representation of axisymmetric surfaces. We use a time-dependent coordinate transformation, which allows us to obtain the shape dynamics in an implicit surface representation. The explicit coordinates of surface points can be calculated independently.

Using our approach, we have identified a mechanochemical shape instability of a sphere that leads to concentration and flow patterns with a polar asymmetry and an accordingly asymmetric oblate shape. For periodic cylinder surfaces, we find contractility-induced instabilities beyond a critical value together with the Plateau–Rayleigh instability for contractility parameters smaller than . Near the point that separates these two regimes, we find that the instability is oscillatory. This could result from a competition between a weak contractile instability and the stability of the cylinder with respect to the Plateau–Rayleigh criterion in this parameter region. On tubular surfaces with stabilized neck radii, we have described the spontaneous formation of directed surface flows. In a reference frame of vanishing average surface flows, this corresponds to a propagating surface constriction that resembles the dynamics of a peristaltic wave (Movie S5).

Our framework provides a basis to explore a large variety of systems that involve the mechanochemical self-organization of deforming active surfaces. Here, we have focused on simple cases, where the active material is embedded in an environment with a homogeneous pressure. The latter exerts a force per area normal to the surface with position-independent magnitude. In many real situations, the surrounding medium has material properties that can give rise to more complex patterns of external forces, including tangential shear forces acting on the surface. Furthermore, it will be interesting to consider active surface material properties with different constitutive relations, such as viscoelastic systems and materials that are anisotropic or chiral. Biological examples of such materials are given by epithelial tissues with planar polarity (1) and by anisotropic cytoskeletal systems (30). Furthermore, the minimal model of a single stress regulator chosen here could be extended to more complex chemical schemes that involve several species and interactions between them. These extended models could represent, for example, the actin dynamics present in the cell cortex and its biochemical regulators (3) or the behavior of sets of morphogens spreading in epithelial tissues (11). Finally, we note that many biological processes, such as the migration of cells or invagination events in tissues (13, 27), lead to shapes that are not axisymmetric. A generalization of our approach to nonaxisymmetric surfaces therefore provides an important challenge.