A minimal physical model captures the shapes of crawling cells

Cell motility in higher organisms (eukaryotes) is crucial to biological functions ranging from wound healing to immune response, and also implicated in diseases such as cancer. For cells crawling on hard surfaces, significant insights into motility have been gained from experiments replicating such motion in vitro. Such experiments show that crawling uses a combination of actin treadmilling (polymerization), which pushes the front of a cell forward, and myosin-induced stress (contractility), which retracts the rear. Here we present a simplified physical model of a crawling cell, consisting of a droplet of active polar fluid with contractility throughout, but treadmilling connected to a thin layer near the supporting wall. The model shows a variety of shapes and/or motility regimes, some closely resembling cases seen experimentally. Our work strongly supports the view that cellular motility exploits autonomous physical mechanisms whose operation does not need continuous regulatory effort.

Cellular self-propulsion in 3D gels or tissues might involve either 'grabbing' of the extracellular matrix [2], or 'swimming' via actomyosin flow in the cell interior [5,6,7]. But, since 3D motility is hard to study experimentally, detailed mechanisms remain largely speculative. Much better understood is a case we address here: cells moving on a 2D solid substrate.
The crawling motion of such cells combines two processes [2]. One is actin polymerization, in which monomers add to one end of an actin filament and leave from the other. Given a net polarization of the actin network, this 'treadmilling' pushes the leading edge of the cell forwards. The second process is the 'walking' of myosin motors along filaments. This creates a contractile [8] stress that, alongside interfacial tension and/or membrane elasticity, opposes cellular elongation. Contractility thus assists retraction of the rear part of a crawling cell (which includes the nucleus). Also implicated in crawling are assemblies of intracellular membrane proteins called focal adhesions [1], which enhance mechanical connectivity between actomoyosin and the solid substrate. These increase the effective wall friction and improve crawling efficiency (just as a child can crawl more easily on a pile carpet than a polished floor). Focal adhesions also signal alterations to the balance between treadmilling and contractility [9,10]. Notably though, significant motility is maintained even in laser-ablated or reconstituted cell fragments that lack most or all of the cell's regulatory machinery [11].
Our goal here is to set up a minimal model of crawling motility, to see how fully the observed behaviour can be understood using physics alone. To achieve this goal, we intentionally rely on a very simplified description of complex biochemical and biophysical processes such as F-actin polymerisation and actomyosin contractility. We also employ only an elementary description of the cell membrane, as well as of its coupling to the interior and exterior of the cell. The advantages of such an approach is that there are arguably as few free parameters as possible within the model, and this allows us to isolate the key physical elements of cell crawling (and/or spreading [12]) at surfaces. In other words, we can directly answer the question posed above: can cellular motility be seen as a physics-based machinery which is coupled to, and controlled by, a complicated underlying biochemical network, but is nevertheless able to function independently? The alternative is that the biochemical detail is fundamental, and no aspects of it can be disregarded, even for a basic understanding of motility. This is an important question, both to understand the basic mechanism of cell motility, and to inform subsequent biomimetic work intending to replicate cellular behaviour with synthetic components. In this respect, our simulation model can be viewed as the in silico analogue of in vitro experiments with cell extracts and cell fragments [11].
The drawback inherent within our approach is that, from the beginning, some of the important biological detail is lost. One can argue that, if experimental observations are reproduced by the minimal framework, this detail may be inessential to the fundamental biophysics of motility, and this by itself is very useful information. Nevertheless, in order to check the picture and to prove the previous statement unambiguously, one should eventually insert back realistic details one by one, and explicitly compare the results obtained with those of the fully simplified, minimal model. We initiate this last approach here as well, for a few cell biological features that our baseline model does not capture well. These include the allowance for a relatively high membrane viscosity, and for wall slip, or equivalently, retrograde flow [13]; see Supplementary Note 1, 2, 3 for details of these checks. In our minimal model, we represent a cell or cell extract as a droplet of active polar fluid (actomyosin) immersed in a second fluid and confined by interfacial tension [7]. (The latter is similar to, but simpler than, the fixed area constraint imposed by the outer membrane of a real cell.) As in real cells, actin treadmilling is confined both to the leading edge of the cell, and within a finite thickness λ from the substrate. Our droplet interacts with a solid wall through a controllable friction, represented by a partial-slip boundary condition. In contrast to most previous computer models for cell motility [7,14,15,16,17,18,19,20,21], ours treats cells as 3D objects. Previous 2D models often incorporated a number of important biological details, such as variable adhesion with the substrate [15], or the dynamics of signalling proteins [18]. A difference with our framework, besides the general philosophy highlighted above, is that these 2D models typically view the cell from above and cannot resolve the flow within the cell, or the presence or absence of slip at the boundary, whereas our model can resolve these features. Importantly, our 3D framework also allows us to directly address the factors influencing the lamellipodium shape, which is the most commonly observed morphology in crawling cells. We are aware of one previous 3D study [22] which recreated the motion of flat cell fragments. Ref. [22] differs from our approach in its aims, as it focuses on a specific type of cell motility, rather than attempting to describe morphological transitions between non-motile and motile cells or between different kinds of motility, as we do here.

Results
Crawling dynamics in quasi-2D. To begin, we consider a quasi-2D geometry describing a thin slice through a crawling droplet (Fig. 1a). At the wall we assume parallel anchoring of actin polarization (represented by black arrows in the figure), and also a no-slip boundary condition on fluid flow as expected when focal adhesions are plentiful (the effect of a small slip is considered below). In these quasi-2D simulations actin polarization (hence treadmilling) is confined both to the leading edge of the cell, and within a finite thickness λ from the substrate. This is to model real cells where actin filaments are present, and polarized, mainly in a cortex towards the exterior [23]. Importantly, no qualitative differences in the results are observed if we choose instead a simpler representation of treadmilling, where this is no longer localised at the leading edge, but still confined within a thickness λ from the substrate, see Supplementary Fig. 1. At the droplet's interface with the surrounding fluid we impose outward normal anchoring of variable strength β. The equations of motion are as previously published [7], and detailed in the Methods Section. Fig. 1a shows two snapshots of such a quasi-2D droplet in steady state. Upon increasing the treadmilling rate w 0 , we observe progression from a droplet that moves along the wall with only modest distortion, to one that extends a thin protrusion at its leading edge (also see Supplementary Movies 1 and 2). In both regimes the droplet moves with a speed close to w 0 . Fig. 1b plots the steady-state wall contact area against w 0 , in the case of no-slip boundary conditions (s = 0, red/plusses curve) and partial slip (s = 0.02, s = 0.04, s = 0.06, green/crosses, blue/asterisks and purple/squares curves respectively). The parameter s controls the amount of slip of the fluid at the solid wall, and is defined such that s = 0 corresponds to no-slip while s = 1 corresponds to full-slip (see Supplementary Note 1 for its mathematical definition). The parameter s is thus inversely related to the density of focal adhesions between the substrate and the cell. In both cases, with or without wall slip, the contact area increases sharply near a critical value w c (s), suggestive of a dynamical transition or sharp crossover between two regimes. This transition becomes sharper as the amount of slip on the substrate augments (green/crosses curve), and discontinuous at critical w c for higher values of s (blue/asterisks and purple/squares curves). The strong hysteresis of the last curve confirms its discontinuous nature. Note that in these quasi-2D simulations, we have set actomyosin contractility to zero. This suggests that actomyosin contraction is not essential to the crawling mechanism; indeed, including ζ leads to almost unchanged results (see Supplementary Fig. 1). Hence the protruded shape observed in crawling cells can indeed be achieved by just a balance of actin treadmilling against interfacial tension at the cell membrane. Finally in Fig. 1c, the crawling speed of our droplet (in steady state) is plotted against the slip parameter s while the treadmilling rate is fixed. As expected, increasing the amount of slip strongly reduces droplet motility via the crawling mechanism. Additionally, we still observe a discontinuous transition at critical slip s c above which no protrusion is observed any longer. This concurs with a previous comparison of crawling with contractility-driven 'cell swimming', in which the wall friction was crudely represented by a depth-integrated drag term [7].
Crawling dynamics in 3D. The quasi-2D simulations in Fig. 1 give intriguing hints to the physics of lamellipodial protrusions, which are the pancake like structures often seen at the leading edge of crawling cells. However, the structure of such cells is fully three-dimensional, and the lamellipodium is only one of several morphologies attained by live cells, whether crawling or stationary, when attached to a solid support. We have therefore performed full 3D simulations of a slight variation of the same minimal model as we used in 2D. With respect to the 2D case, in the 3D case the effect of actomyosin contractility is much stronger, so we present results only with it switched on from the beginning. Unlike actin treadmilling, actomyosin contractility is known to show highest activity deep inside the cell, close to the cell nucleus, rather than at the cell periphery [23]. To simplify the model in 3D, actin polarization (and thus contractility) is therefore now assumed to be uniform throughout the droplet. Additionally, actin polarization is no longer confined to the leading edge of the cell, but we retain a spatially varying treadmilling rate, w 0 → w(z) = w 0 exp(−z/λ), to confine actin treadmilling close to the substrate. The use of this simplified model, where the active droplet is homogeneously polarised, is motivated by our finding above that the localisation of actin polymerisation is not essential to observe the morphological transition between non-motile and motile cells in Fig. 1 (see Supplementary Fig. 1).
As we shall see, our minimal 3D model supports the formation of various dynamical cell shapes, some with surprisingly direct experimental counterparts. Fig. 2 shows a selection of these steady state morphologies, together with a direct comparison with experimental morphologies taken from the literature. Fig. 2a-c shows a droplet subject to strong normal anchoring of the po-larity (large β). This boundary condition, which might model the promoted nucleation of outwardly polarized actin by Arp2/3 complexes close to the outer membrane [2], enforces radial polarization with an 'aster' defect at the centre [24]. Treadmilling then causes the droplet to spread, but not to move. In steady state it adopts a symmetric 'fried egg' whose shape resembles that of an inactive platelet [1] or a stationary keratocyte [25]. (We speculate that adding boundary pinning of the cell perimeter, and/or radial modulation of the treadmilling strength w 0 , might create shapes more closely resembling mesenchymal stem cells [26]).
Figs. 2d-f, 2g-i, 2j-l show steady-state cell shapes for (rightward) crawling cells; the morphology is selected by a combination of treadmilling rate w 0 , contractility −ζ, and anchoring strength β. Fig. 2d-f (see also Supplementary Movie 3) arises when all three effects compete; this closely resembles a motile keratocyte (see e.g. [2,27]). The pronounced lamellipodium at the leading edge stems from a balance between contractility-induced splay of the polarity field [7,8], and its soft anchoring at the leading edge. The prediction of this familiar dynamical structure from such minimal physical ingredients is a major finding of this paper. (In contrast the overhang at the cell rear is a detail, and can be eliminated if actin polarization or treadmilling rate w(x, y, z) is set to zero at the back of the cell, as done e.g. in Fig. 1a.) The cell shape in Fig. 2g-i arises when the anchoring strength β is decreased further: contractility-induced splay now dominates, creating an indented structure possibly reminiscent of the 'phagocytic cup' which forms when a cell needs to engulf a solid particle [28]. (Such cups are typically not planar but have been reported for crawling as well as free cells [29,30].) Finally, Fig. 2j-l shows a droplet where treadmilling is dominant: both the anchoring and contractility are much weaker than in 2d-f or 2g-i. Importantly, without those contributions the cell shape no longer shows the keratocyte-like lamellipodium morphology. Instead, the leading-edge forms a pointed protrusion, qualitatively resembling a 'pseudopod', as often reported in ameoboid cell motion both on surfaces and in 3D geometries [31,32].
Phase diagram. The results shown in Fig. 2 were obtained by varying the treadmilling parameter w 0 , the contractile activity −ζ, and also the anchoring β which controls the orientation of the polarisation at the cell boundary. Although a systematic study of the phase diagram in this multidimensional parameter space is extremely demanding computationally, and hence outside the scope of the current work, we observe (see Fig. 3) that similar morphologies can also be obtained by keeping the value of w 0 constant, varying β, and adjusting ζ to be as large as possible avoiding at the same time droplet breakup (which occurs when ζ is too large; see also Ref. [33] where a distinct but related study of contractility-driven breakup 2D active droplets is presented).
For very small values of β, a finger-like structure (the "pseudopodium" in Fig. 2) is dominant, while fan-shaped geometries (such as the "lamellipodium" in Fig. 2) emerge when β is larger, due to the higher splay distortion of the actin filaments arising at the leading edge of the cell. For very strong anchoring, the lamellipodium morphology becomes unstable, and the cell attains a fried-egg shape, which is non-motile. Experiments (see e.g. Ref. [27]) have studied how cell shape is affected by adhesion to the substrate. In our model an increase in the adhesion would affect both the treadmilling and the anchoring, decreasing the first and increasing the second. Interestingly, our simulations suggest that lamellipodia are found for intermediate values of the anchoring strength and are disfavoured when treadmilling decreases, qualitatively confirming what has been experimentally observed by Barnhart et al. [27].
Evolution of polarisation and velocity fields. Besides the concentration field whose structure is shown in Fig. 2, our model also tracks the evolution of the polarization and intracellular flow fields, P and u respectively.
In Fig. 4 we show the polarization field close to the substrate (left column) and the flow field (right column) corresponding to the 3D steady state morphologies reported in Fig. 2. The polarization field shows a transition between the aster corresponding to the 'fried-egg' state ( Fig. 4a) and polarized uniform (Fig. 4d) and splayed (Fig. 4b,c) pattern; these are qualitatively in broad agreement with the experimental data for keratocytes presented in [25], where it was established that motility requires symmetry breaking of the actin fibre orientation field, from an aster like conformation as in our "fried-egg" morphology ( Fig. 4a) to an asymmetric conformation in a lamellipodium, where the barbed ends of actin fibres are approximately normal to the cell membrane at the leading edge (as in our simulations in Fig. 4c).
Coming to the flow field, we first note that u (plotted in Fig. 4) should be thought of as an intracellular solvent flow inside the cell, which has not been measured often in experiments (an exception is Ref. [34] which is discussed in the Supplementary Note 2). This intracellular flow is quite distinct from the flow of F-actin, which is usually tracked in cell biology experiments, and which is commonly referred to as "retrograde flow" [13]. F-actin flows inwards within the cell, from the leading edge to the bulk. The velocity field for the "fried-egg" cell is characterized by a quadrupolar velocity profile (in a fully 3D view), more intense close to the surface of the droplet where the curvature is higher than near the substrate (see Fig. 4a). The flow in the vicinity of the membrane is not unlike that measured in the case of an activitydriven motile tumour cell "swimming" within a 3D matrigel [5], and even closer correspondence is found between that experiment and our contractilitydriven droplets (which move in the absence of treadmilling, see Fig. 5). The symmetric structure of the flow field in the "fried-egg" morphology explains its non-motile character. The flow field for the keratocyte, lamellipodiumlike shape, and for the phagocytic cup ( Fig. 4b-c) is instead unidirectional (along the direction of motion) in the laboratory frame of reference, and gently follows the fan-out structure of the surface cell profile, being more intense close to both cytoplasmic projections. A similar structure is seen for the "pseudopodium" (Fig. 4d), with the difference that the flow is now almost entirely localised around the finger-like protrusion. Interestingly the highest magnitude of the fluid flow for the last three cases is twice larger than the one for the fried egg cell.
The flow fields of these structures in are consistent with previous computer simulations [35], and with the flow fields mapped experimentally in Ref. [34] with decreased contractility compared to the wild-type (see also Supplementary Fig. 2). Experiments suggest that contractility can propel intracellular flow significantly, and 2D simulations with variable ζ in our model confirm this qualitatively (see Supplementary Fig. 3 and Supplementary Note 2 for further discussion).
Oscillatory dynamics and medium-dependent viscosity. Not all parameters in our model lead to steady state shapes; we also find spontaneous oscillations in shape and motility (see Supplementary Fig. 4, Supplementary Note 3, and Supplementary Movie 4). These or similar oscillations have long been reported experimentally [36] and it is remarkable that they can emerge in principle from a model that is devoid of all biochemical feedbacks.
All our simulations address the case where the active droplet is surrounded by a fluid of similar viscosity to its interior. This is an adequate model of wall-bound cells moving through an aqueous medium as is normally the case in vitro and in vivo, where cells are surrounded by at least a thin film of fluid. Accounting for the higher viscosity of the membrane also leads to similar results; we have performed selective numerical checks of this ( Supplementary Fig. 5).
Crawling in the absence of treadmilling. Motility of the crawling droplets in Fig. 2d-f, 2g-i, 2j-l is powered by treadmilling, although contractility was found to play a crucial role in determining cell shape. What happens if treadmilling is absent altogether? It was shown in [7] that when suspended in bulk fluid a contractile droplet can 'swim' by undergoing spontaneous splay. This creates an internal toroidal fluid flow which entrains the external fluid and leads to self-propulsion. Our simulations show that the same physical mechanism is still available for a fluid-immersed droplet attached to a wall (Fig. 4). The flow associated with such contractility-driven active droplets is similar to that observed experimentally in cells 'swimming' in 3D matrigel (Fig. 5).

Discussion
In conclusion, we have presented a minimal physical model of cell crawling, based on an active fluid droplet, attached to a solid wall and immersed in a second fluid. Key ingredients in the model are actin treadmilling (confined to a thin layer of depth λ close to the substrate), actomyosin contractility (throughout the droplet), full anchoring of the actin polarity parallel to the wall, and its partial anchoring normal to the fluid-fluid interface. Hydrodynamic fluid flow is properly treated throughout. We showed that on increasing the treadmilling rate, such a droplet undergoes a nonequilibrium morphological transition to form a protrusion layer whose thickness is set by λ. Variation of the treadmilling, contractility and anchoring parameters lead in 3D to a striking series of cell morphologies, resembling structures seen in real eukaryotic cells including both lamellipodia and pseudopods. We also find regimes where oscillatory modulations of shape and motility arise. Our results show that surprisingly many of the observed features of motile and spreading cells (though of course not all of them) can emerge from relatively simple physical ingredients. This suggests that some aspects of cellular motility might fruitfully be thought of as driven by a physics-based 'motility engine' whose function, although controlled by the cell's complex biochemical feedback networks, does not directly depend upon these for its basic operational principles. It would be of interest to extend the current model to study cell motility inside complex media, such as the extracellular matrix, or within tissues [37].

Methods
We briefly outline here the hydrodynamic model used in this work. We consider a droplet of cellular matter sitting on a solid substrate as a simple representation of a crawling cell. The droplet phase is represented by a scalar field φ(r, t) > 0. Outside the droplet, we have a passive isotropic fluid with φ(r, t) = 0 to represent the wet environment around the cell and the substrate. Secondly we define a polarization field P(r, t) to represent the actin filaments which are polarized and abundant mainly at the cortex of the cell. More specifically, P is the coarse grained average of all the orientations of the actin filaments: P = p where p is a unit tangent oriented from the end at which the actin filaments depolymerize to the opposite end at which they polymerize. Finally, we also have the average velocity field v(r, t) of both the cellular matter and the solvent.
Underlying free energy. Following the phenomenological approach used in Refs. [7] and [38], we first introduce a free-energy functional that controls the equilibrium physics of a polar liquid-crystalline passive droplet The first two terms of Eq. (1) control bulk and interfacial properties of the droplet. In particular, k > 0 determines the droplet interfacial tension and a > 0 is a phenomenological constant. The bulk term is chosen to stabilise a droplet phase (φ φ 0 ) immersed in a second fluid (φ = 0) and φ cr = φ 0 /2. The parameter t(r) controls a second order transition from isotropic (|P| = 0) to polar phase (|P| > 0) inside the droplet as the value of t changes from 0 to 1. Note the for all our 3D work, and for selected cases in 2D (see Supplementary Fig. 1), we have taken t(r) = 1 everywhere for simplicity; we discuss some implications of this choice below. The functional form of t(r) for the remaining quasi-2D simulations is: where z is the vertical distance from the substrate, y is the coordinate parallel to the substrate and y back is the position of the cell rear. Basically, the value of t is equal to 1 when z is less than some thickness λ from the substrate and at least a distance Λ from the cell rear. This defines the polarized (hence active) region at the leading edge of the droplet where actin treadmilling takes place (see black arrows in Fig. 1A). Incidentally, λ also sets the lengthscale for the thickness of the protrusion layer. Finally, elsewhere in the droplet t is equal to 0 which gives rise to an isotropic phase (and hence no actin treadmilling).
In summary, there are three equilibrium phases which can be obtained by minimizing the free energy F [φ, P] in a state of uniform φ and P, setting δF δφ = 0 and δF δP = 0. These three phases are: the isotropic passive fluid, external to the droplet (φ = 0 and P = 0); the polarized region of the droplet, corresponding to the cortical actin-myosin network (φ = φ eq1 and P = P eq ); and the isotropic region of the droplet, corresponding to cell nucleus and other cellular materials which do not play a role in cell motility (φ = φ eq2 and P = 0).
The other three terms in the free energy describe bulk and elastic properties of the polar field, where α is a positive (for stability) constant and κ is the elastic constant (in the one-elastic constant approximation [39]) that describes the energetic cost due to elastic deformation in a liquid-crystalline environment. The surface tension γ will depend on a, k and on the elastic constant κ. Across the droplet interface, the values of the concentration and the polarization fields vary smoothly from φ = φ eq and P = P eq to φ = 0 and P = 0. In general the width of the diffuse interface also depends on the same parameters. In our simulations, the interfacial width is typically much smaller than the size of the droplet. Finally the last term βP · ∇φ takes into account the anchoring of P to the droplet surface (so that for β > 0 the polar field P preferentially points perpendicularly outwards at the droplet perimeter while for β < 0 it preferentially points perpedicularly inwards).
Equations of motion. The velocity field v(r, t) obeys the continuity and the Navier-Stokes equations, which in the incompressible limit are Here P is the isotropic pressure and σ active + σ passive is the total stress, given by the sum of active and passive contributions. The first of them is that can be derived by summing the contributions from each force dipole and coarse graining [38]. The constant ζ is the activity parameter which is negative for contractile particles and positive for extensile particles. Its magnitude |ζ| is proportional to the strength of the force dipole. For the case of actomyosin networks, ζ is negative signifying the tendency to contract along the direction of actin filaments. The strength of actomyosin contraction is thus given by its magnitude |ζ|. As we have seen from Fig. 2, actomyosin contractility is responsible for the formation of a variety of complex 3D structures in crawling cells. Also note that in our calculations above, we neglect the contributions of higher order terms (such as ∂ α P β ) in the active stress, although these are in principle allowed by symmetry [40]. This assumption is appropriate for actomyosin solutions, which are made up by elongated fibres of uniform width whose degree of asymmetry between the particle's head and tail is supposed to be small.
We note that actomyosin contraction shows highest activity at the cell rear [23] rather than at the leading edge. This is in contrast to actin treadmilling which shows highest activity at the leading edge of the cell. Thus to study the global effect of actomyosin contractility in 3D, we make a further simplification by assuming actin polarization to be homogenous and uniform inside the droplet, which amounts to setting t(r) = 1 throughout. Since actin polarization is now no longer localised at the leading edge, we instead localise the treadmilling parameter w to be large near the substrate to produce a protrusion as in the quasi-2D systems. While these approximations are useful to make progress in 3D, as the numerics is more demanding there, we note that the protrusion transition found in Fig. 1 is still preserved, both in 3D and in 2D when similar approximations are made ( Supplementary Fig. 1).
The passive stress consists of a viscous term, written as σ viscous αβ = η(∂ α v β + ∂ β v α ), where η is the shear viscosity (the greek indices denote cartesian components), summed with an elastic stress σ elastic , due to elastic deformations in the liquid crystal phase, and a surface tension term σ interf ace . The elastic stress is given by similar to that used in nematic liquid crystal dynamics. Here the constant ξ depends on the geometry of the active particles: ξ > 0 describes rod-like particles, ξ < 0 oblate ones. We have set it positive, as in [7]. In addition ξ controls whether particles are flow aligning in shear flow (|ξ| > 1), hence creating a stable response, or flow tumbling (|ξ| < 1), which gives an unstable response. Here we assume ξ > 1. Note that the elastic stress also depends on the molecular field h, defined as δF/δP. The surface tension term takes into account interfacial contribution between the passive and the active phase, and is given by This term is borrowed from binary fluids [41], with f being the free energy density.
The dynamics of the concentration of active material φ(r, t) is governed by a convective-diffusion equation where M is a thermodynamic mobility parameter, related to the diffusion constant D M a (see Eq. (1)). δF/δφ is the chemical potential of the system and w 0 is the self-advection parameter. This w 0 can be interpreted as the speed, relative to the bulk fluid velocity v(r, t), at which each filament is self propelled along its own tangent [7]. Thus w 0 is also proportional to the rate of actin treadmilling. In our 3D minimal model, the treadmilling rate is a spatially varying function which is large close to the substrate: w 0 → w(r) = w 0 exp(−z/λ) where λ is the lengthscale for protrusion thickness.
The last ingredient to complete the description of the physics of the system is an evolution equation for the polarization field P(r, t). This is given by where Γ is the rotational viscosity and D = (W + W T )/2 and Ω = (W − W T )/2 are the symmetric and the anti-symmetric part of the velocity gradient tensor W αβ = ∂ α v β . Note that this equation is the one usually employed in polar liquid crystal theory [42].
To solve these equations we use a 3D hybrid lattice Boltzmann algorithm already successfully tested in other systems, such as binary fluids [43], liquid crystals [44] and active matter [7]. It consists in solving Eq.(8) and Eq.(9) via a finite difference predictor-corrector algorithm while the integration of the Navier-Stokes equation is taken care of by a standard Lattice Boltzmann approach. This numerical approach is necessary due to the couplings between velocity field and polarization field in Eq.(9) and velocity field and concentration in Eq.(8).

Supplementary Note 1: Boundary conditions and the effect of wall slip
To simulate a cell crawling on a flat solid surface, we need to set up an appropriate treatment of droplet wetting at the substrate. To this end, we first initialized our system with a half-spherical droplet on the plane, with concentration φ = φ eq inside and φ = 0 outside. The polarization field inside the droplet is initially uniform and points along the y-direction, while it is zero outside the droplet. The droplet is therefore initially polarised; the system then evolves via the equations of motion described in the previous section.
At the substrate the polarization P has a parallel and strong anchoring, i.e. P z (z = 0) = 0 and ∂Px ∂z z=0 = ∂Py ∂z z=0 = 0, while the anchoring of the polarization at the surface of the droplet is provided by the term βP · ∇φ in the free energy (see the main text). With β > 0, the polarization points outwards and is perpendicular to the surface of the droplet, while with β < 0 it points inwards (remaining perpendicular to the surface). For the actomyosin concentration φ at the wall we adopt neutral wetting which means that a · ∇φ| z=0 = 0 and a · ∇(∇ 2 φ)| z=0 = 0, where a is a unit vector normal to the substrate. These conditions together ensure a · ∇µ| z=0 = 0, so that the concentration gradient is parallel to the wall and there is no flux across it. Importantly, focal adhesions, which prevent the filaments from slipping, are taken into account by using non-slip boundary conditions at the substrate (z = 0). This means that v(z = 0) = 0. Focal adhesions play a fundamental role in cells sliding on a surface. Their absence can cause the actin filaments to slip backwards as they polymerise against the cell membrane. On the other hand, experimentalists can reduce focal adhesions by coating the glass slide with appropriate enzymes before putting the cell (such as a keratocyte) on top [45]. In order to emulate the reduced presence or lack of focal adhesions at the surface, we have also simulated a system in which there is a mixture of non-slip and full-slip boundary conditions (see e.g. Fig. 1 in the main text, and Supplementary Fig. 1). This is done by following the approach proposed in Ref. [46] in which a parameter s is introduced to model the "fraction of the slip" at the wall. More specifically, a partial slip condition is defined to where i = x, y and s is the slip parameter (s → 0 corresponds to non-slip condition while s → 1 corresponds to full-slip condition). The "density" of focal adhesion between the cell and the substrate might be quantified via the term 1/s -if there are fewer focal adhesions the filaments are more likely to slip. Fig. 1C in the main text shows a plot of the crawling speed as a function of the slip parameter for a 2D simulation. Since the slip velocity at the surface will be in the opposite direction to that of the actin polymerization, an increase in the amount of slip causes a decrease in the crawling speed. This is consistent with the fact that the actin polymerization becomes less effective at cell protrusion with fewer focal adhesions [45,2].
Supplementary Note 2: Flow field in the frame of reference of the cell The flow fields in Fig. 4 of the main text are measured in the reference frame of the laboratory. Within the frame of reference of the cell, the flow is mainly backwards (Supplementary Fig. 2), in agreement with the simulations of Ref. [35] and with the experimental results in Ref. [34] for the case of a cell treated with blebbistatin. Contractility experimentally was found to be able to reverse the sense of this flow. In qualitative agreeement with this finding, we observe that contractility leads to a larger flow (see the intracellular flow in the 2D simulations in Supplementary Fig. 3, where simulations with wall slip and with ζ = 0 or ζ = 0 are directly compared). Finally, it is possible to make semiquantitative contact with experiments measuring retrograde actin flow by tracking, instead of u, u + wP. In reality, P is non-zero only close to the surface, as in our 2D simulations in Fig. 1 of the main text (see also Supplementary Movies 1 and 2). Supplementary Fig. 3 shows the intracellular flow field for 2D simulations, in the frame of reference of the cell, which equals u − V cell (V cell wP is the cell velocity). In the same frame of reference, the F-actin flow would be equal to u + wP − V cell u (V cell wP is the cell velocity). Close to the surface (z ∼ 0), where the experiments typically track retrograde flow, the F-actin flow velocity is therefore approximately u(z = 0), hence there is no retrograde flow unless there is wall slip (as in the case of Supplementary Fig. 3). This should be expected, because strong focal adhesion leads to no slip, strong friction, hence no retrograde actin flow.

Supplementary Note 3: Oscillatory dynamics and medium-dependent viscosity
Our minimal physical model supports the formation of several 3D cell morphologies, some of them remarkably close to their experimental counterparts (see Fig. 2 in the main text). We recall here that in our 3D simulations, unlike in typical 2D simulations (e.g. in Fig. 1 of the main text), we have simplified the actomyosin contractility to be uniform throughout the droplet (t(r) = 1). However, in our 3D simulations actin polarization is confined to a thin layer close to the substrate (unlike the 2D case in Fig. 1, where actin polarization and treadmilling are only effective at the leading edge of the cell). This is achieved by assuming a spatially varying treadmilling whose functional form is w(z) = w 0 exp(−z/λ), where z is the Cartesian coordinate perpendicular to the substrate and λ is related to the thickness of the leading edge. 2D simulations using this prescription lead to a similar transition as observed in Fig. 1 of the main text (see Supplementary Fig. 1.) While the structures shown in Fig. 2 and Fig. 3 are the most common ones, as noted in the main text, these are not the only ones supported by our model. In particular, for suitable choices of parameters the steady state is no longer a stable crawling morphology as described in the main text, but can also become oscillatory.
In Supplementary Fig. 4 (a-h) we show the time evolution of one such oscillatory motion, where the active droplet glides on a flat surface while its leading edge continuously switches between two different shapes. A typical crawling mechanism is initially observed (a) while, later on, the leading edge of the droplet undergoes a deformation characterized by a depression (b) that splits the lamellipodium into two "feet" (c and d). These then separate and later on gradually merge (e); the droplet then crawls again with a single lamellipodium and eventually the cycle repeats (f , g, h). This dynamics is qualitatively reminiscent of the oscillatory crawling motion of keratocytes on substrates with high friction [45]. Similar sequences of shape changes are also observed when cells develop a phagocytic cup before engulfing foreign bodies such as solid particles [28]. Our results therefore show that similar shape changes can be due to cell elasticity, actin treadmilling and motor activity alone.
Finally, we come back to the fact that our simulations in the simplest version address a droplet moving within a fluid of its same viscosity. In reality, the two fluids will not have the same viscosity (the cell usually moves within a thin film of an aqueous medium), and the membrane viscosity will be significantly larger than that of the cytosol. To simulate this, we consider a case in which the viscosity is intermediate in the cell (the cytoplasm viscosity is ∼ 10 cP [47]), high within the membrane, and low outside. This more complex situation leads again to our realistic cellular morphologies, an example with a lamellipodium-like protrusion is shown in Supplementary Fig. 5.

Supplementary Note 4: Mapping to physical units
In this section we briefly review how to map simulation to physical units in our work (see also Ref. [7] for more details); this mapping is recapped in Supplementary Table 1. Parameters have been chosen to be of the same order of magnitude as the ones quoted in [45,35]. In our case, realistic parameter ranges can be achieved by fixing the length-scale, time-scale and force-scale as: L = 0.3 µm, T = 10 ms, and F = 1000 nN respectively (in simulation units these scales are all equal to one exactly). This correspondence then determines the mapping to physical units of all quantities.
Note that, in line with standard practice in Lattice Boltzmann simulations, the fluid mass density ρ is much larger than the actual mass density of a real solvent (water) [48]. This is acceptable so long as inertial effects remain sufficiently small, and this procedure speeds up the computations by several orders of magnitude. In this case the choice of the force/density scale above keeps the Reynolds number below 0.1 which is small enough for inertial effects to be negligible. Figure 1: Quasi-2D crawling dynamics. (a) shows a quasi-2D actomyosin droplet crawling along the treadmilling direction (to the right in the picture). The black arrows here represent actin polarization which points in the same direction as actin treadmilling. Additionally, actin polarization is localised at the leading edge of the droplet to model a crawling eukaryotic cell (also see Supplementary Movies 1 and 2). If the treadmilling rate w 0 is large enough, a thin protrusion forms, whose steady-state size (area per unit transverse length) is plotted versus w 0 in (b) (red/plus curve), for the case of no-slip boundary condition between the droplet and the substrate. From (b), we can identify the critical treadmilling rate w c above which a thin layer of protrusion is observed in a crawling cell. This transition becomes sharper as the amount of slip at the substrate increases (green/cross curve), and discontinuous/first order for higher values of s (blue/asterisk and purple/square curves). The hysteresis loop corroborates its discontinuous character. (c) shows a plot of the crawling speed as a function of the slip parameter s (while the treadmilling rate is fixed at w 0 = 0.002). This controls the discontinuity in fluid velocity at the solid wall, and is defined such that s = 0 corresponds to no-slip while s = 1 corresponds to full-slip. An increase of s causes a decrease of the velocity of the cell. This is because the slip velocity at the surface is in the opposite direction to that of the actin polymerization. Again, we observe a discontinuity at s = s c marking a dynamical phase transition from a round shaped crawling cell (s > s c ) to a protruded shape (s < s c ). Figure 2: Spreading and crawling cell morphologies. In this figure spreading and crawling cell morphologies observed in past experiments (left column) and obtained in the current 3D simulations (middle and right column, giving respectively top and side view) are shown. Frame a-c shows a static cell with a fried egg shape; the left column panel was originally published in Ref. [25]. Other frames show moving cells with protrusions resembling respectively, a lamellipodium (d-f; the left column panel was originally published in Ref. [27]), a phagocytotic cup (g-i; left column panel taken from Ref. [30] with permission), and a pseudopod (j-l; left column panel taken from Ref. [32] with permission). Parameters for the 3D modelling (in simulation units; see Supplementary Table 1)   , lamellipodium (c) and fried egg crown (d)). Notice that w 0 is constant; while ζ varies (both ζ and β affect droplet elasticity, increasing ζ for small β leads to further instabilities such as droplet breakup which are not relevant for the current work).  Figure 5: Flow field of a contractile active droplet. The figure shows the flow field (in the laboratory frame of reference) for a contractilitydriven motile active droplet (left from Ref. [7]; for simplicity in the bulk), and for a cell moving in 3D inside matrigel (from Ref. [5], with permission). Figure 6: 3D dynamics in absence of treadmilling. Snapshots corresponding to the motion of a contractile cell on a substrate, in the absence of treadmilling (parameters were w = 0, ζ = −0.005, β = 0). The active droplet first elongates perpendicular to the wall due to the active stress; the polarization then splays in 3D and this leads to internal flow and a motile concave shape, resembling half of a free-swimming 3D droplet as investigated in [7].