Viscous vertex model for active epithelial tissues

TL;DR

The paper develops a rotation-invariant, viscous extension of the vertex model that incorporates both junctional viscosity along cell–cell interfaces and bulk viscosity between vertices and cell centers. It regularizes the zero-friction limit with a Lagrange-multiplier formalism, enabling well-posed simulations and flexible boundary conditions, and introduces a slab-shear protocol to extract a coarse-grained tissue viscosity ${\eta}_{\rm tissue}$ from microscopic parameters. Analytically, ${\eta}_{\rm tissue}^{(ST)} = {\eta_s}/{(4\sqrt{3})} + {\eta_b}/{(2\sqrt{3})}$ for a hexagonal cell, and numerically, long-time viscosity under sustained shear scales with ${\eta_s}$ and ${\eta_b}$ as cell rearrangements occur. When active polar or nematic stresses are added, increasing cellular viscosity elongates cells, reduces defect density, and promotes coherent flows, providing a direct link between cell-resolved dissipative physics and continuum active-nematic descriptions in free-floating tissues and organoids. The framework thus offers a practical bridge between discrete cell models and continuum rheology, with broad relevance to epithelial mechanics and morphogenesis.

Abstract

We present a rotationally invariant viscous vertex model that accounts for both cortical and bulk dissipations of cells. The vanishing substrate-friction limit is enforced via Lagrange multipliers, which also provide a route to strongly constrained boundary conditions such as fixed boundaries and prescribed tractions. Building on this formulation, we introduce a slab-shear rheology protocol to extract an effective, coarse-grained tissue shear viscosity. Under polar or nematic activity, viscosity regulates the formation of elongated, spatially correlated cell-shape textures and stabilizes well-defined topological defects. Because the model remains well-posed at zero substrate friction, it is naturally suited to describing free-floating epithelia and organoids.

Figures (6)

• Figure 1: Schematic of the active, viscous vertex model. (a) In vertex models, a cell sheet is mimicked by a tilting of polygons. (b) Schematic of the two kinds of cellular viscosities considered in our model, including the cell--cell interfacial viscosity ($\eta_{ij}^{(s)}$) and the cell cytoplasm viscosity ($\eta_J^{(b)}$). (c) Schematic of the cell--cell interfacial tension of the interface $ij$, including three contributions: (1) a constant tension $\Lambda_{ij}^{(0)}$; (2) a viscous tension $\eta_{ij}^{(s)}\dot{l}_{ij}$; (3) tension fluctuations $\delta \Lambda_{ij}$. (d) Schematic of the cell--bulk line tension of the segment $iJ$, including three contributions: (1) a constant tension $\Lambda_{J}^{(0)}$; (2) a viscous tension $\eta_{J}^{(b)}\dot{l}_{i,J}$; (3) tension fluctuations $\delta \Lambda_{i,J}$. (e) Schematic of the polar active traction force $\bm{F}_J^{(\rm active)} = T_0 \bm{p}_J$ with $T_0$ being the active traction magnitude and $\bm{p}_J$ the cell polarization vector. (f) Schematic of the apolar active stress $\bm{\sigma}^{\rm (active)}_J = -\beta \bm{Q}_J$ with $\beta$ being the cellular activity parameter and $\bm{Q}_J$ the cell shape anisotropy tensor.
• Figure 2: Numerical simulation of an active cell sheet in a square domain with periodic boundary conditions. Here, we consider cell-shape-dependent active stresses; see Sec. \ref{['sec:activity_nematic']}. (a) The condition number of the friction-viscosity coefficient matrix $\bm{C}$ as a function of the friction $\gamma$, for a system consisting of $N_c = 10^{3}$ cells. The dashed blue line represents the condition number of the extended viscosity coefficient matrix $\bm{C}_{\rm ext}$ (see Eq. \ref{['eq:C_extended']}). Data analysis gives a scaling law, ${\rm cond}(\bm{C}) \sim \gamma^{-1}$. (b-d) The case without friction ($\gamma = 0$), using the protocol proposed in Sec. \ref{['sec:zero-friction']}. (b) The condition number of the extended viscosity coefficient matrix $\bm{C}_{\rm ext}$ (see Eq. \ref{['eq:C_extended']}) as a function of the number of cells $N_c$. Data analysis gives a scaling law, ${\rm cond}(\bm{C}_{\rm ext}) \sim N_c^{\alpha}$ with $\alpha \approx 2.6$. (c) A typical cell morphology. The black lines indicate cell orientations; the red (resp. blue) symbols represent $+1/2$ (resp. $-1/2$) topological defects, extracted using the scheme proposed in Ref. Lin2023. (d) A typical flow field, where the black arrows represent the velocity vectors and the color code refers to the velocity magnitude. (e, f) Floppy modes of vertices' motions in the case of no friction ($\gamma = 0$) and no bulk viscosity ($\eta_b = 0$). (e) Illustration of a typical floppy mode in a cell layer consisting of $N_c = 100$ cells, obtained by numerical calculation of Eq. \ref{['eq:floppy_mode']}. The arrows represent the vertices' velocity vectors. (f) The fraction of floppy modes $f$ as a function of the number of cells $N_c$. Parameters: $\eta_s = 10$, $T_0 = 0$, and $\beta = 0.4$; $\eta_b = 10$ for (a-d).
• Figure 3: Measuring the coarse-grained tissue viscosity. Here, we do not consider cell activities (i.e., $T_0 = 0$ and $\beta = 0$). (a) Schematic of shearing a cell sheet within a slab geometry. The slab is assumed to be periodic along the $x$-axis and is of width $L_y$ in the $y$-axis. Cell vertices adhered to the bottom border are fixed, i.e. $\bm{v} = \bm{0}$, while cell vertices adhered to the top border move at a specified velocity $\bm{v} = v_0\bm{e}_x$. Thus, the global tissue shear strain rate is $\dot{\gamma}_{xy} = v_0 / L_y$. (b) Theoretical analysis of a single hexagonal cell under a constant shear strain rate $\dot{\gamma}_{xy}$. (c) The short-time tissue viscosity $\eta_{\rm tissue}^{(\rm ST)}$ and the long-time tissue viscosity $\eta_{\rm tissue}^{\rm (LT)}$ as functions of cellular viscosity $\eta = \eta_s = \eta_b$. Comparison with the analytical result (Eq. \ref{['eq:eta_tissue_analytical']}). (d-g) Comparisons of the cell deformation strain field (d, f) and the cell flow field (e, g) for the case with friction (d, e) and the case without friction (f, g), averaged over $N_t = 100$ frames. In (d, f), the color code refers to the magnitude of the cell deformation tensor $\bm{\varepsilon}_{\rm cell}^{(\rm dev)}$ (the deviatoric part) and lines represent orientations. In (e, g), the color code indicates the magnitude of velocity, and the arrows represent the velocity vectors.
• Figure 4: Numerical simulation of pulling a cell within a passive cell sheet with a constant pulling velocity $v_0$. Here, we do not consider cellular activities and assume $\eta_s = \eta_b = \eta$. Comparison of (a, b) the pure frictional case ($\eta = 0$) and (c, d) the strong viscous case ($\eta = 100$). (a, c) Evolution of cell morphologies. Here, we mark the cell under pulling in green. (b, d) Evolution of the isotropic stress within cells. (e) Comparison of the average cell elongation (of neighboring cells of the cell under pulling) for the pure frictional case and the strong viscous case. (f) Comparison of the average isotropic stress (of neighboring cells of the cell under pulling) for the pure frictional case and the strong viscous case. (g) Comparison of the pulling force for the pure frictional case and the strong viscous case. Parameters: $T_0 = 0$, $\beta = 0$, $v_0 = 0.01$, and $\gamma = 1$.
• Figure 5: Numerical simulation of viscous, polar, active tissue flows in a circular domain. (a-c) The cell morphology (left) and the cell velocity field (right) at various levels of cell viscosity $\eta_b = \eta_b = \eta$: (a) $\eta = 1$; (b) $\eta = 10$; (c) $\eta = 100$. The black lines represent cell orientation and the red (resp. blue) symbols indicate the locations and orientations of $+1/2$ (resp. $-1/2$) topological defects, extracted using the scheme proposed in Ref. Lin2023. The arrows represent cell velocity vectors, and the color code refers to the cell velocity magnitude. (d) The number of topological defects and the rotational order parameter as functions of the cellular viscosity $\eta$, averaged over $n = 5$ independent simulations. Parameters: $T_0 = 0.05$, $D_r = 0.05$, $\mu_{\rm LA} = 0.05$, and $\mu_{\rm CIL} = 1$.