Poisson-bracket formulation of the dynamics of fluids of deformable particles

Many extended systems, such as biological tissue lecuit2007cell , foams durian1995foam ; cohen2013flow , emulsions mattsson2009soft ; vlassopoulos2014tunable , and colloidal suspensions mattsson2009soft can be described as collections of deformable particles. A variety of mesoscopic models have been developed to examine the role of particle shape on the structure and rheology of these soft materials.

Cellular Potts models graner1992simulation ; kabla2012collective and Vertex and Voronoi models honda2004three ; hufnagel2007mechanism ; farhadifar2007influence have been successfully used to describe dry foams and confluent layers of biological tissue, where cells completely cover the plane with no gaps, with extensions to three dimensions hannezo2014theory ; murisic2015discrete . These models describe cells in confluent tissues as tightly packed irregular polygons covering the plane and predict a jamming-unjamming transition tuned by a target cell shape that captures the interplay of cortex contractility and cell-cell adhesion, with the mean cell shape serving as a metric for tissue fluidity bi2014energy ; bi2015density ; bi2016motility . Vertex and Voronoi models do not, however, have a natural extension to situations where the cell packing fraction is below one, although gaps between cells have been incorporated in recent work teomy2018confluent ; kim2020embryonic . In contrast, both particle deformability and density variations can be incorporated in multi-phase field models and in models of deformable particles boromand2019role , which have been used to examine solid-liquid transitions as a functions of both particle shape and density.

Less well developed are continuum descriptions of the rheology of materials where the constituents can change their shape. An important example is the classic work by Doi and Ohta that describes the dynamics of the interface between two immiscible fluids under shear, incorporating formation, rupture and deformation of droplets doi1991dynamics . Continuum mechanics of confluent tissue have been constructed phenomenologically and employed to connect structure and mechanics in Drosophila development sagner2012establishment ; popovic2017 . Ishihara and collaborators formulated a continuum model that couples cell shape to mechanical deformations at the tissue scale ishihara2017cells . Their work, however, only captures simultaneous cell anisotropy and alignment of elongated cell shapes, without distinguishing between a tissue where cell shapes are on average isotropic and one where cells are on average anisotropic, but not aligned, as observed in simulations of Vertex/Voronoi models bi2016motility ; yang2017correlating . It is in fact the single-cell anisotropy that provides an order parameter for cell jamming in Vertex and Voronoi models bi2014energy ; bi2015density ; bi2016motility , where fluid states of elongated cells are obtained without nematic order of elongated cells. The importance of this distinction in a continuum theory of tissue mechanics was highlighted recently in work by one of us and collaborators czajkowski2018hydrodynamics .

In this paper we adopt the Poisson-bracket formulation forster1974microscopic to obtain continuum equations for a fluid of deformable particles in two dimensions. This method has the advantage of providing a systematic derivation of the reversible part of the hydrodynamic equations once the continuum fields have been identified. Our approach is inspired by work by Stark and Lubensky Stark2003 ; Stark_2005 who used the Poisson-bracket approach to derive the hydrodynamics of a nematic liquid crystal. As in liquid crystals, we identify both a continuum scalar field that quantifies fluctuations of individual cell shape and a cell shape tensor field that captures both cell elongation and alignment. An important difference is that, while in passive liquid crystals molecular shape fluctuations decay on fast (non-hydrodynamic) time scales, in a tissue cell shape is the order parameter for the rigidity transitions, hence cell-shape fluctuations are long-lived near the transition and must be incorporated in a hydrodynamic model. The equations derived here provide a continuum model for collections of interacting deformable “particles” and can be adapted to describe both confluent and non-confluent systems.

The paper is organized as follows. In Section II we provide the microscopic definition of the continuum fields used in the hydrodynamic model. In Section III we briefly summarize the Poisson Bracket (PB) method and the calculation of the various PBs (with details given in Appendix C), and discuss the reactive and dissipative contributions to the coarse-grained dynamics. The final continuum equations are displayed in Section IV. In Section V we discuss the form of the continuum equations for the specific case of a cellular tissue, and conclude with a brief discussion in Section VI. Details of the derivation of the PBs and of the mean-field free-energy of the Vertex model are given in Appendices.

We consider a fluid whose constituents are $N$ extended particles of arbitrary shape. The contour of each particle, referred to below as a ‘cell’, is described by a polygonal shape joining $n$ vertices located at $\mathbf{r}_{\mu}^{\alpha}$, where $\mu=1,2,\cdots,n$ labels the vertices and $\alpha=1,2,\cdots,N$ labels the cells, as shown in Fig. 1. Each cell has a total mass $m_{c}$, which we assume equally distributed among the $n$ vertices. The shape of each cell is described by a shape tensor defined as

where $\Delta\mathbf{r}^{\alpha\mu}=\mathbf{r}^{\alpha\mu}-\mathbf{r}_{\alpha}$, with $\mathbf{r}_{\alpha}=\frac{1}{n}\sum_{\mu}\mathbf{r}^{\alpha\mu}$, and Latin indices $i,j$ denote components.

We define microscopic mass, momentum and cell shape density fields as

with $m=m_{c}/n$. Coarse grained quantities are then defined as $\rho(\mathbf{r},t)=\left[\hat{\rho}(\mathbf{r},t)\right]_{c}$, $\mathbf{g}(\mathbf{r},t)=\left[\hat{\mathbf{g}}(\mathbf{r},t)\right]_{c}$ and $G_{ij}(\mathbf{r},t)=\left[\hat{G}_{ij}(\mathbf{r},t)\right]_{c}$ and correspond to macroscopic continuum fields describing the system on length scales large compared to both the size of the particles and their mean separation. Note that since the microscopic cell-shape tensor $\hat{\mathbf{G}}^{\alpha}$ has dimensions of length squared, the density of cellular shape tensor $G_{ij}$ is dimensionless. As we will see below, the trace of the shape tensor density provides a measure of the density of cell perimeter, while its traceless part, $\tilde{G}_{ij}=G_{ij}-\frac{1}{2}\delta_{ij}\text{Tr}[\mathbf{G}]$, captures both cell anisotropy and local alignment of elongated cells.

The cellular shape tensor can be written in terms of its eigenvalues as

where $\lambda_{1}^{\alpha}>\lambda_{2}^{\alpha}$ and $\hat{\bm{\nu}}^{\alpha}$ is the eigenvector of the largest eigenvalue. For regular $n-$sided polygons, the shape tensor is diagonal with $\lambda_{1}^{\alpha}=\lambda_{2}^{\alpha}$. In this case the cell area $A_{\alpha}^{(n)}$ and perimeter $P_{\alpha}^{(n)}$ can be expressed in terms of the invariants of the tensor $\mathbf{G}^{\alpha}$ as

Single cell anisotropy is measured by $\Delta_{\alpha}=\lambda_{1}^{\alpha}-\lambda_{2}^{\alpha}$ which vanishes for regular polygons. To quantify single-cell elongation independently of alignment of elongated cells, we follow czajkowski2018hydrodynamics , albeit with a slightly different definition of the shape tensor, and introduce an anisotropy density field defined as

and the associated coarse grained field $M(\mathbf{r},t)=[\hat{M}(\mathbf{r},t)]_{c}$. Work on Vertex/Voronoi models of confluent biological tissue, as well as multiphase fields models, has demonstrated the correlation between tissue fluidity and anisotropy of single cell shape, as quantified here by $M$. In Vertex models, this anisotropy provides an order parameter for the solid-liquid transition bi2015density ; bi2016motility .

In the following, we construct hydrodynamic equations for a fluid of deformable particles that couple structural changes encoded in cell shape and alignment of elongated cells to flow. The dynamics of the fluid on scales large compared to the cell size and mean cell separation is described in terms of a few continuum fields: the mass density $\rho$, the momentum density $\mathbf{g}$, the single-cell anisotropy density $M$ and the cell-shape tensor density $G_{ij}$.

Here we briefly summarize the Poisson-Bracket (PB) formalism. Consider a system whose microscopic dynamics is determined by canonically conjugate positions $\mathbf{r}^{\alpha}$ and momenta $\mathbf{p}^{\alpha}$. We describe the dynamics in terms of a few microscopic density fields $\hat{\Psi}^{a}(\mathbf{r},t;\{\mathbf{r}^{\alpha}\},\{\mathbf{p}^{\alpha}\})$, for $a=1,2,\cdots$. These fields are chosen to be be either hydrodynamic fields associated with conserved quantities, broken symmetry fields, or quasi-hydrodynamic fields that decay on times scales large compared to microscopic ones. In the specific case of interest here $\{\hat{\Psi}^{a}\}=(\hat{\rho},\mathbf{\hat{g}},\hat{G}_{ij},\hat{M})$. The dynamics of the corresponding coarse-grained fields $\Psi^{a}(\mathbf{r},t)=[\hat{\Psi}^{a}(\mathbf{r},t;\{\mathbf{r}^{\alpha}\},\{\mathbf{p}^{\alpha}\})]_{c}$ is governed by the equations

where $V^{a}$ and $D^{a}$ represent the non-dissipative and dissipative parts of the dynamics, respectively. The reactive term $V^{a}$ is given by

where $\mathcal{F}[\{\Psi^{a}\}]$ is a phenomenological free energy,

and

Finally, the dissipative term in the kinetic equation is controlled by all the neglected microscopic degrees of freedom and can be written as

The dissipation tensor $\Gamma^{ab}$ is in general a functional of the $\{\Psi^{a}\}$ and their gradients. It is a phenomenological quantity controlled by the requirement that $\partial_{t}\Psi^{a}$ can only couple to driving forces $\frac{\delta\mathcal{F}}{\delta\Psi^{b}(\mathbf{r})}$ that have different sign under time reversal, to guarantee that such terms describe dissipation. In equilibrium it is a symmetric tensor and must obey Onsager’s principle de1951thermodynamics .

Putting it all together, we now write the final form of the equations we have obtained. It is convenient to write

where $D_{ij}$ is the rate of strain tensor given in Eq. (32) and $\omega_{ij}$ is the vorticity,

The set of continuum equations for our fluid of deformable cells is then given by

where

are the convective derivative, and the convective and corotational derivative.

The equation for the shape tensor $G_{ij}$ contains couplings to flow vorticity and strain rate which control the tendency of extended and deformable particles to rotate with flow and align with streamlines. The shape tensor $G_{ij}$ plays a role similar to that of the conformation tensor in a polymer suspension beris1994thermodynamics . In fact, if we ignore the additional anisotropy density field $M$, the equations derived here for a fluid of deformable particles have the same structure as a one-fluid model of viscoelastic polymer solutions milner1993dynamical . Unlike in models of polymer suspensions, however, the coefficient of the coupling to strain rate, known in that context as the slip parameter beris1994thermodynamics , is found to be simply equal to $1$ in our PB formulation.

It is also convenient to separate the dynamics of $G_{ij}$ in that of its trace and deviatoric part. The corresponding equations are given by

where $[...]^{ST}$ denotes the symmetrized and traceless part of the tensor. These equations need to be completed by an expression for the free energy $\mathcal{F}_{V}$ in terms of the shape tensor and anisotropy density field. Such an expression of course depends on the system of interest. In the next section we consider the specific case of a model of dense biological tissues.

Confluent biological tissue, where cells are tightly packed, with no intervening gaps, have been modeled extensively using Vertex or Voronoi models that describe cells as irregular polygons tesselating the plane honda2004three ; hufnagel2007mechanism ; farhadifar2007influence . The behavior of the tissue is controlled by an energy that describes the tendency of each cell to adjust its area $A_{a}$ and perimeter $P_{a}$ to target values $A_{0}$ and $P_{0}$, given by

with $\kappa_{A}$ and $\kappa_{P}$ stiffness parameters. The first term arises from tissue incompressibility in three dimensions and the second captures the interplay of cell-cell adhesion and cortical contractility. By scaling lengths with $\sqrt{A_{0}}$ and energies with $\kappa_{A}A_{0}^{2}$, the scaled energy of each cell is given by

with $p_{0}=P_{0}/\sqrt{A_{0}}$ the target shape parameter and $r=\kappa_{P}/(\kappa_{A}A_{0})$.

Numerical studies of this energy have identified a rigidity transition at a critical value $p_{0}^{*}$ of the target shape parameter between a rigid, solid-like state for $p_{0}<p_{0}^{*}$ to a fluid state for $p_{0}>p_{0}^{*}$. Single-cell anisotropy as quantified by the mean cell-shape index $q=\langle P_{\alpha}/\sqrt{A_{\alpha}}\rangle$, with the brackets denoting an average over all cells, provides an order parameter for the transition. Czajkowski et al. czajkowski2018hydrodynamics derived a mean-field model of this rigidity transition, albeit using a different definition of the cell shape tensor $G_{ij}^{\alpha}$ as compared to the one used here. The derivation carried out with our definition is outlined in Appendix  D. The result is a quartic Landau-type free energy density $f_{M}$ where the the cell shape anisotropy density $M$ plays the role of an order parameter, given by

where $\alpha(p_{0})$ vanishes at $p_{0}=p_{0}^{*}$ and $\beta>0$. The definition of the shape tensor of individual cells only affects the precise values of these parameters that also depend on the reference polygonal shape, but does not change the form of the free energy density nor the value of $p_{0}^{*}$. The free energy given in Eq. (49) is obtained by assuming small deformations from regular polygons and constant cell perimeter. It predicts a mean-field transition at $\alpha=0$ from a state where cells are isotropic ($M=0$) $\alpha<0$ or $p_{0}<p_{0}^{*}$ (the solid state) to a state where cells are anisotropic ($M=\sqrt{-\alpha/\beta}$) for $\alpha<0$ or $p_{0}>p_{0}^{*}$ (the liquid state).

This work suggests a phenomenological free energy for a confluent tissue that captures both fluctuations in the cell anisotropy density $M$ that quantifies the liquid-solid transition and the shape tensor density $\tilde{G}_{ij}$ that quantifies alignment of elongated cell as

We do not include terms of order $Tr[\tilde{\mathbf{G}}^{2}]^{2}$ as we do not expect any nematic order of cellular shapes in the absence of externally applied or actively generated internal stresses. Also, we have assumed constant cell perimeter, corresponding to $R=Tr[\mathbf{G}]=$constant. In general, the various parameters in $\mathcal{F}_{c}$ will depend on $R$.

It is important to stress that $\tilde{G}_{ij}$ and $M$ are not independent. The traceless tensor $\tilde{G}_{ij}$ can be written as

which defines the director field $\mathbf{n}(\mathbf{r},t)$ associated with alignment of elongated cells and the magnitude $S$ of orientational order. Cell alignment can only occur if cells are elongated ($M\not=0$), hence $S_{G}(M)$ must vanish when $M=0$. We assume $S_{G}=MS$, where $S$ plays the role of a nematic order parameter for orientational order of elongated cells. Clearly, $S$ is defined only in states where $M$ is finite.

Cell sheets commonly interact with a frictional substrate that eliminates momentum conservation. Frictional drag with the substrate generally exceeds inertial forces, and the Navier-Stokes equation for the momentum is replaced by a Stokes equation quantifying force balance on each fluid element. Within this overdamped limit, and considering a minimal form for the various dissipative kinetic coefficients, the tissue dynamics is governed by the following equations

where $\xi$ is the frictional drag and

It is useful to consider a simplified form of the equations obtained by retaining only lowest order terms in fields and gradients. In this case the Stokes equation and the equations for the shape fields can be written in the explicit form

where $D=\Gamma_{M}K$, $D_{G}=\Gamma K_{G}$, $r=\chi/\gamma_{G}+\alpha/\gamma_{2}$ and

The single-cell anisotropy field $M$ here plays the role of tissue fluidity. The first term on the RHS of Eq. (59) captures the fact that shear deformations, coupled to local cell alignment, can increase cell anisotropy, driving fluidification. The second term describes relaxation to the ground state controlled by the tissue free energy, with a cost for spatial variations in local fluidity controlled by the stiffness $D$. The reactive terms in Eq. (60) describe flow alignment of elongated cell shape. The term proportional to $r$ describes changes of cell shape tensor due to dissipative processes, such as topological rearrangements, at a rate proportional to the tissue fluidity $M$. The last term in Eq. (60) describes the stiffness against deformations of local cell alignment.

Finally, in a confluent tissue the cell number density $n=\rho/m_{c}$ is slaved to the mean cell area $\langle A_{\alpha}\rangle$ with $n=1/\langle A_{\alpha}\rangle$. For cells that are only slightly deformed from regular polygons, $\langle A_{\alpha}\rangle\simeq\sqrt{\det[\mathbf{G}]}\approx\text{Tr}[\mathbf{G}]$, where we have used Eq. (66). The density equation, Eq (52), can therefore equivalently be written as an equation for the cell area or for $|G|\equiv\det[\mathbf{G}]$, given by

Using the Poisson bracket formalism, we have derived hydrodynamic equations for a fluid of deformable particles in two dimensions. Shape fluctuations are described by two continuum fields: (i) a coarse-grained scalar field that captures single-particle anisotropy, and (ii) a shape tensor field that quantifies both particle elongation and nematic alignment of elongated particles.

We have specifically applied the model to sheets of dense biological tissue, where single-cell anisotropy was recently identified as the order parameter for a solid-liquid transition driven by the interplay of cortex contractility and cell-cell adhesion bi2015density ; bi2016motility . In other words, in confluent tissue single-cell anisotropy is effectively an experimentally accessible measure of the rheological properties of the tissue, with isotropic cell shapes identifying the solid or jammed state and anisotropic shapes corresponding to a liquid. Previous work has examined the dynamics of a coarse-grained cell shape tensor and its coupling to mechanical stresses ishihara2017cells . This work did not, however, distinguish between a tissue of elongated, but isotropically oriented cells and one were the cells and elongated and also aligned in a state with nematic liquid crystalline order. The new ingredient of our work is to distinguish the dynamics of tissue fluidity, as quantified by the single-cell anisotropy field, from that of cell alignment, and examine the interplay between flow, which can be either externally applied or induced by internal active processes, fluidity and nematic order of cell shapes. Our equations hence provide a starting point for quantifying the rheology of biological tissue. Future extension needed to develop a complete framework of tissue rheology include the coupling to the dynamics of polarized cell motility and the inclusion of structural rearrangements arising from cell division and death.

Finally, the equations developed here provide a general hydrodynamic model for any fluid of deformable particles, capable of accounting for both the dynamics of shape deformations and density changes.