{'fraction_non_alphanumeric': 0.03955025150409311, 'fraction_numerical': 0.014301213137390275, 'mean_word_length': 4.346809632624363, 'pattern_counts': {'":': 0, '<': 6, '<?xml version=': 0, '>': 6, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Vertex model instabilities for tissues subject to cellular activity or applied stresses

The vertex model is widely used to describe the dynamics of epithelial tissues, because of its simplicity and versatility and the direct inclusion of biophysical parameters. Here, it is shown that quite generally, when cells modify their equilibrium perimeter due to their activity, or the tissue is subject to external stresses, the tissue becomes unstable with deformations that couple pure-shear or deviatoric modes, with rotation and expansion modes. For short times, these instabilities deform cells increasing their ellipticity while, at longer times, cells become non-convex, indicating that the vertex model ceases to be a valid description for tissues under these conditions. The agreement between the analytic calculations performed for a regular hexagonal tissue and the simulations of disordered tissues is excellent due to the homogenization of the tissue at long wavelengths.

The vertex model is widely used to describe the dynamics of epithelial tissues, because of its simplicity and versatility and the direct inclusion of biophysical parameters. Here, it is shown that quite generally, when cells modify their equilibrium perimeter due to their activity, or the tissue is subject to external stresses, the tissue becomes unstable with deformations that couple pure-shear or deviatoric modes, with rotation and expansion modes. For short times, these instabilities deform cells increasing their ellipticity while, at longer times, cells become non-convex, indicating that the vertex model ceases to be a valid description for tissues under these conditions. The agreement between the analytic calculations performed for a regular hexagonal tissue and the simulations of disordered tissues is excellent due to the homogenization of the tissue at long wavelengths.

I. INTRODUCTION
The vertex model, initially proposed to describe foams and soap bubbles [1,2], has been extended to describe epithelial tissues [3][4][5][6] with large success. Applications include the study of cell division [7], tissue elongation [8] and epithelial packing in wing disk and ventral furrow formation in Drosophila [9][10][11], tube formation [12,13], and the rigidity transition in active tissues [14]. Approximating each c cell as a polygon, an energy functional is built that penalizes the deviations of the actual cell areas (A c ) and perimeters (P c ) from preferred values (A 0c and P 0c , respectively). In the most generic form, the energy functional is with l i j the length of the cell edge shared by vertices i and j. K A is the area elastic modulus, which describes the three dimensional incompressibility of the layer and the resistance to height fluctuations; K P is the perimeter elastic modulus related to the actin-myosin contractility; J is the adhesion energy per unit length and represents a constant line tension. Although it is possible to absorb the last term into the second one by redefining P 0c , we opt to keep all terms, such that the different constants retain a direct interpretation. Through this work we consider A 0c and P 0c given by the initial geometry of each cell. Hence, the model only has three free parameters.
In its usual form, the degrees of freedom of the model are the positions of the vertices r i , which evolve variationally as where γ is a mobility that we will absorb in K A , K P and J, which now have units of relaxation rates times different powers of length. Active stresses are continuously induced by cell divisions, extrusions and rearrangements between neighboring cells [15]. Also, stresses are generated by cell growth [16] and contractions [17]; processes that can be easily included in the vertex model as changes in the equilibrium cell parameters.
In Refs. [5,10], the vertex model was used to obtain the phase diagram of the ground state (the most relaxed network configuration) of a proliferating tissue, initially made of a regular hexagonal packing. They found a phase transition induced by cell division in the parameter space One phase corresponds to a single ground state, with regular hexagonal packing geometry, while the other phase corresponds to a network with many soft deformation modes, where the hexagonal packing looses stability. Here, we develop a general framework to study the stability of tissues subject to cell activity and externally applied stresses. Neither cell division nor cell rearrangements are considered. This is the case of some experiments [18,19] and previous analytical calculations [5,20,21]. Also, topological events are non-linear and, therefore, they are not relevant to describe the emergence of the instabilities. We show that for a large region of the parameter space, if in large portions of the tissue the cells modify their activity or it is subject to external stresses, the whole tissue becomes unstable in the form of long-wavelength deformations that couple pure-shear or deviatoric modes, with rotation and expansion modes. These instabilities differ from those that take place in passive foams [22,23], because they are triggered by the cellular activity.
The organization of the paper is as follows. In Sec. II we present the general analysis of the instabilities that appear in a confluent tissue, focusing in the case of cellular activity. The analytical method for regular tissues and the simulations for irregular ones are described and compared. Section III considers the case of tissues subject to external pre-stresses. In Sec. IV we discuss the case of general anisotropic pre-stresses, which need a more detailed analysis. Our conclusions and a discussion are presented in Sec. V. Finally, the Appendices give technical details.

II. TISSUE UNDER CELL ACTIVITY
For the vertex model, the elastic coefficients K A and K P are assumed to be positive, and penalize deviations from the reference areas and perimeters, while there is no restriction on the sign of J, as has been discussed in the literature [6,24]. As a first case, where analytical results can be obtained, we consider a regular tissue composed of N identical hexagonal cells of side a, for which A 0c = 3 √ 3a 2 /2 and P 0c = 6a, for arXiv:2010.12482v1 [cond-mat.soft] 23 Oct 2020 2 all cells c. Cell activity can generate stresses that tend to deform the tissue. For example, sudden changes in the actomyosin activity in the cell border can be modeled as a modification of the equilibrium perimeters, P 0c → (1 + λ P ) P 0c (with λ P > 0 for expansions and λ P < 0 for contractions). Similarly, a change in the actomyosin activity in the medioapical side of the cells imply changes in the equilibrium cell areas, As a first case, we consider homogeneous modifications of the tissue (uniform λ P and λ A ), modeling large portions of the tissue that change as in Ref. [19], and we investigate the stability and rigidity of this tissue, allowing it to fluctuate. The vertex positions are now given by (I + εU) r 1, and U a general 2 × 2 matrix of components u ik , characterizing the fluctuations. Computing contributions up to O(ε 2 ), the energy of the tissue may be written as , where the superscripts represent the order of each term in the expansion, and E A , E P and E J are the contributions proportional to K A , K P and J, respectively. The full expressions are given in Appendix A 1.

The stress tensor is
. It has a zeroth order contri- pλ P δ ik that represents the total stress, with passive and active contributions, needed to maintain the deformed configuration. Here, we defined the energy scaleÊ = NK A A 2 0 /2 and the dimensionless parameters p = K P /(a 2 K A ) and j = J/(a 3 K A ), which are the ratios between the characteristic time of the surface elasticity and the ones related to the perimeter and adhesion elasticity, respectively.
For general fluctuations, U can be expanded in Fourier modes. When computing the total energy of the tissue, the linear terms in ε cancel by spatial integration, leaving only the reference energy and the quadratic terms in the fluctuations. In physical terms, the linear contribution is eliminated by the application of a uniform external stress σ ik by other tissues that act as a frame, imposing rigid boundary conditions. Furthermore, in the limit of small wavevectors k, the dominant contribution comes from the case of homogeneous U, plus small corrections proportional to k 2 , which we neglect henceforth. Hence, to analyze the stability of the tissue under long wavelength fluctuations, we have to determine whether the quadratic form for homogeneous U is positive definite. Expressing U as a linear combination of four basic deformation modes, In the case where the deformation is due to cell activity, the µ-matrix is diagonal with where we used the expressions of Appendix A 1. The deformation modes U 1 and U 2 are both shears, although in different directions. Consequently, their eigenvalues, which are associated to the shear modulus, are equal. Negative values of the diagonal terms signal the development of an instability of the corresponding mode, in a single cell description. For example, large positive values of λ A (cell expansion), would give rise to unstable rotation and expansion modes, while for large negative values of λ A (cell compression), the deviatoric and pure shear modes become unstable. At a tissue level, however, due to the confluent property, pure modes are not allowed. Indeed, consider for example the Fourier mode where the new vertex positions are given by x = x + ε sin(2πx/L) cos(2πy/L) and y = y − ε sin(2πy/L) cos(2πx/L), shown in Fig. 1a. Depending on the position, some cells experience deviatoric deformations (in yellow), while others rotate (in red). Similarly, for the Fourier mode x = x + ε cos(2πx/L) sin(2πy/L) and y = y + ε cos(2πy/L) sin(2πx/L), shown in Fig. 1b, pure shear modes (in green) coexist with expansion modes (in blue). Simple uniaxial deformations with a sinusoidal amplitudes also couple the deviatoric and expansion modes. Complementary to the long wavelength fluctuations, it is possible that the boundaries between neighboring cells move inside a supercell (analogous to optical phonons in solids) as shown in Figs. 1c and d. Again, different modes coexist. The confluent property with the periodic boundary conditions frustrate the emergence of pure deformation modes. The use of fixed boundary conditions leads to the same frustration.
This unavoidable coexistence of modes implies that even though a deformation mode may seem to be unstable at the cell level, the total energy of the tissue should be computed as the sum of the different contributions that, at the end, may result to be positive definite. A detailed study of the stability of a tissue that considers the coexistence of modes is given in Section IV. We provide here a qualitative argument to obtain the stability limit from the behavior of individual cells. As the deviatoric and pure shear modes share the same value in the µ-matrix, the total energy of the tissue fluctuations shown in Figs. 1a and c are equal, with a prefactor equal to µ 11 + µ 33 = j/3−4pλ P . An instability is hence predicted to develop for λ P > j/(12p). Notably, when λ A = 0, the instability is predicted to take place when the shear modulus (i.e. µ 11 or µ 22 ) vanishes, as was observed in Ref. [5]. However, when the target area has changed (λ A = 0), the vanishing of the shear modulus does not signal the development of unstable modes.
To validate the predictions in actual situations, we simulate both regular and irregular tissues. Regular hexagonal tissues are made of N = 3000 cells arranged in a box of size L x = 50 √ 3a and L y = 90a with periodic boundary conditions. In order to avoid artificial effects due to the lattice perfections, a Gaussian noise is added to all the vertex po- nsion], , although ues which ative valf an instascription. n), would while for deviatoric ent propsider for positions y 0 = y ent property with the periodic boundary conditions frustrate the emergence of pure deformation modes. The use of fixed boundary conditions leads to the same frustration. This unavoidable coexistence of modes implies that even though a deformation mode may seem to be unstable at the cell level, the total energy of the tissue should be computed as the sum of the different contributions that, at the end, may result to be positive definite. A detailed study of the stability of a tissue where the cells are deformed by activity is given in the Supplementary Material, we provide here a more qualitative argument to obtain the stability limit from the behavior of individual cells. As the deviatoric and pure shear modes share the same value in the µ-matrix, the total energy of the tissue fluctuations shown in Figs. 1a and c are equal, with a prefactor equal to µ 11 + µ 33 = j/3 4p f . An instability is hence pre- x = x + ε sin(2πx/L) cos(2πy/L) and y = y − ε sin(2πy/L) cos(2πx/L), and b) x = x + ε cos(2πx/L) sin(2πy/L) and y = y + ε cos(2πy/L) sin(2πx/L). c) and d) fluctuations where the boundaries between neighbor cells move inside a supercell. Cells with well defined deformation modes are colored: yellow for deviatoric, green for pure shear, red for rotation, and blue for expansion. For simplicity, square cells are used in the presentation. sitions in both directions, with standard deviation 0.1a. Irregular tissues are built as Voronoi cells, where the positions of N = 3000 center points are generated by a Montecarlo simulation of hard disks in a box of equal size as for the regular tissue. The diameter of the disks govern the degree of dispersion of the cells. We consider an area fraction φ = 0.71, below the freezing transition, to obtain a reproducible disordered tessellation with moderate dispersion in cell sizes. The irregular tissues are made of polygons of different sizes and number of sides, implying variance in the equilibrium areas and perimeters, A 0c and P 0c . The deviatoric and pure shear modes manifest in the elongation of cells, which we characterize by the flattening parameter c = (a − b)/(a + b), computed for each cell in terms of its principal semiaxis a and b, calculated as the square root of the eigenvalues of the texture ma- where the sum is over the n c vertices conforming the cell, with positions r i , and r c is the center of the cell. Simulations are performed solving numerically the equations of motion (2), which are worked out in the Appendix B [Eqs. (B2), (B8), (B9), and (B10)]. The differential equations are integrated using the Euler integration method, for various values of K P and J, fixing units such that K A = 1 and a = 1. The time step was fixed to dt = 0.005 and we study the system up to t = 0.5.
The change of the standard deviation of the flattening parameter after few time steps for fixed positive perimeter change λ P = 1/2, considering λ A = ±1/2, displays an important increase precisely where the instability is predicted (Figs. 2a and b). The chosen values of λ A,P are consistent in while positive values indicate that non-convex polygons appear. The agreement with the analytical prediction is excellent, both when regular and irregular tissues are simulated. A comparison between regular and irregular tissues is presented in the Appendix D, showing that the instability takes place for the same parameters and the values of the observables agree. Importantly, the line at which the shear modulus vanishesobtained when neglecting the coupling of modes-fails to predict the instability for all tissues (Figs. 2, 5, and 6). For the cases shown in Figs. 1b and d, the energy for the tissue has a prefactor that becomes negative when pλ P > 3/2 + 2p + j/12, requiring an extremely large increase of the equilibrium perimeter, except if j is negative. Consequently, these modes are hardly seen and are hidden by other more unstable modes.
For cells of equal equilibrium area and complete contraction of the perimeter (λ P = −1), the transition line in Refs. [5,10] is reproduced. An important difference with their work is the use of a fixed size box in simulations, generating at long times non-convex polygons instead of soft networks.

III. TISSUE UNDER PRE-STRESS
In addition to cellular activity, the tissue can be subject to a pre-stress generated by the action of neighboring cells or tissues, fixed boundary conditions, an actomyosin network, or the drag by another expanding tissue located in an adjacent layer, causing it to get pre-deformed. To model a pre-stressed tissue, we perform an affine transformation by changing the vertices positions as r   i . As for the cell activity, we consider homogeneous deformations of the tissue (uniform Λ) and perturbations U in the small wavevector limit, and we analyze first the different deformation modes independently, without dealing with their coupling. For an hexagonal cell, it is found that E A =Ê det(Λ) 2 tr(U) 2 + 2det(Λ) (det(Λ) − 1) det(U) . The expressions for E (2) P and E (2) J are more involved but numerically it is found that they are always positive definite for all pre-deformations, when K P and J are positive (see Appendix A 2 for the full expressions). We conclude, then, that negative J could give rise to instabilities for any pre-strain. The case of E A requires more analysis. From the expression for E (2) A , it is found that fluctuations with det(U) = 0 are always stable.
A is diagonal with el- Note that whenever det(Λ) = 0, either µ A11,A22 or µ A33 are negative, giving rise to possible unstable modes. When det(Λ) > 1 (for example, under a pre-expansion), µ A11,A22 are negative and the deviatoric and pure shear modes may be unstable. Also, when 0 < det(Λ) < 1 (for example, under a compression predeformation), µ A33 is negative and the rotation mode may be unstable. To fully determine the stability, we must consider the perimeter and edge contributions to the energy, as well as the mode couplings.
For isotropic pre-strain Λ = (1 + h)I (h > 0 for expansions and −1 < h < 0 for compressions), the complete µ-matrix is diagonal, with The stability of the relevant global mode is, therefore, described by µ 11 + µ 33 = (1 + h)(4hp + j/3), which can become negative for a wide range of parameters when the tissue is under compression. Simulations are performed, using the methods described in Section II, for an isotropic compression of 50%. Figure 3-left shows an excellent agreement with the analytical calculations that predict the instability line at j = 6p. Again, the instability manifests in an increase of the eccentricity and, at longer times, the appearance of non-convex polygons.  I  II  III  I  II  III FIG. 3. Tissue instabilities obtained in simulations of N = 3000 irregular cells in tissues under 50% isotropic contraction (left), and under 60% horizontal contraction plus 40% vertical expansion (right). Same representation as in Fig. 2.

IV. ANISOTROPIC PRE-STRESSES
Finally, in vivo or in vitro tissues are in general subject to anisotropic external deformations [7][8][9], causing the µmatrix to be non-diagonal. The relevant global modes are obtained as follows. For an extended tissue, the fluctuation is expanded in Fourier modes: r = r + ∑ k a k e ik·r . From the Jacobian of this transformation, the local deformation matrix is computed as u αβ (x, y) = ik α a kβ e ik·r . Expanding it as U(x, y) = ∑ 4 i=1 v i (x, y)U i , a local energy density is obtained, e(x, y) = (Ê/L 2 ) ∑ 4 i, j=1 µ i j v i (x, y)v j (x, y). Finally, the total energy of the tissue is where we used that v i (x, y) are linear combinations of the Fourier coefficients a k and that the Fourier modes decouple if the tissue is homogeneous on the large scale. The matrix e αβ is a 2 × 2 matrix with real coefficients.
where we used that the µ-matrix is symmetric. The stability of the tissue, considering the confluent and periodic conditions, is then obtained from the eigenvalues of the e-matrix, which depend only on the directionk of the wavevector. If at least one eigenvalue is negative, the tissue develop long wavelength instabilities. When the µ-matrix is diagonal, and using that µ 11 = µ 22 , it is found that the eigenvalues of e αβ do not depend on θ and they are given by 1 4 (µ 11 + µ 33 ) and 1 4 (µ 11 + µ 44 ), which corroborates the simple analysis for the coupling of modes described in Section II. Anisotropic pre-deformations generate non-diagonal µmatrices, for which some examples are given in the Appendix E. Figure 3-right presents the comparison between simulations and the prediction of the instability using the eigenvalues of the associated e-matrix for a tissue under 60% horizontal contraction plus 40% vertical expansion. The agreement is again excellent when the non-convexity proxy is used. The flattening parameter does not signal the instability because, for this case there is no manifestation in the change of ellipticity as a result of the coupling of all modes. Finally, Fig. 4 shows the results for a tissue that is subject to a pure deviatoric stress or to a pure shear stress.

V. DISCUSSION
Our analysis shows that stressed tissues described by the two-dimensional vertex model present instabilities in which the cells deform to increase their ellipticity, to later become non-convex. These stresses can be generated by the cellular activity when the actin ring on the perimeter of the cells changes its size or they can be external, when the tissue is pre-stressed. In any of these cases the tissue is unstable for a wide range of the model parameters.
The presence of the predicted instabilities is a stringent test of the vertex model to describe biological tissues, which under many conditions are subject to internal and external stresses. For example, in developing tissues, processes like invaginations, cell extrusion and division generate stresses. Uniaxial pulling can be generated by other tissues [15] or driven experimentally [19,21,25]. Also, biochemical signals can alter in large regions the activity of the tissue [19]. These and other configurations, with different external stresses, should be investigated to verify if the predicted instabilities take place and if they can act as seeds to instabilities in developing tissues. In the mechanobiological approach, forces and instabilities launch the tissue transformations during development that are necessary to generate structures and organs [26,27]. If the vertex or similar models correctly describe the tissue dynamics, internal or external stresses can trigger the instabilities described in this letter, which can initiate tissue transformation processes.
In this letter we restricted the analysis to two-dimensional planar dynamics. Further studies are needed to analyze how the deformation modes couple with motion in the third dimension when the planar restriction is removed. For example, buckling instabilities generating wrinkles, could relax stresses instead of generating non-convex polygons.  For the analytic calculations, we consider a regular tissue composed of N identical regular hexagonal cells of side a, for which the preferred cell area and perimeter for all cells are A 0c = 3 √ 3a 2 /2 and P 0c = 6a, respectively.

Tissue under cell activity
Cell activity is included as homogeneous modifications of the equilibrium perimeters, P 0c → (1 + λ P ) P 0c and equilibrium areas A 0c → (1 + λ A ) A 0c , with λ P , λ A > 0 for expansions and λ P , λ A < 0 for contractions.
We define A Then, when considering an activity modulated by λ A , the term of the energy proportional to K A is given by Hence, the zeroth, first, and second order terms of E A are We define P Then, when considering an activity modulated by λ P , the term of the energy proportional to K P is given by The zeroth, first, and second order terms of E P are therefore given by Finally, the adhesion contribution to the energy is where P c (1) is given in Eq. (A6), As a result, the zeroth, first, and second order terms of E J are given by Eqs. (A1) and (A6) can be obtained using Mathematica.

Tissue under stress
Now, we study the same energy contributions, but when the tissue is subject to a homogeneous strain, such that all the vertices change their position as r  i , where Λ is a 2 × 2 matrix that gives account of the pre-deformation.
In a similar way as in the previous section we can define A (1) c and P (1) c , representing the area and perimeter of the cell c, that was initially a regular hexagon with area A 0c and perimeter P 0c , which is now subject to a given strain characterized by the matrix Λ. Then, we define A (2) c and P (2) c as the values when we allow fluctuations, modulated by the matrix U, in the system.
The expressions for P (1) c and P (2) c are more complicated to write in terms of the matrices Λ and U. In general terms, considering that the six vertices of the hexagon have positions r i , we obtain: with where we use r i+1,1 = r i+1 − r i , assuming the vertices ordered clockwise. The terms M (1) c and M (2) c are given by Now, following a similar procedure as in the previous section we can compute all the energy contributions. The contribution proportional to K A is where we obtain that the zeroth, first, and second order terms of E A are given by Similarly, for the term proportional to K P , and the zeroth, first, and second order terms of E P are given by Finally, the zeroth, first, and second order terms of E J are c . (A33)

Appendix B: Equations of motion
With periodic boundary conditions, Eq. (1) from the main text can be written as The equations of motion for the vertex are obtained using Eq. (2) of the main text, which can be written as Assuming a polygon of N vertices, we calculate its area using the triangularization method with respect to the vertex v 1 , where we used that the tissue is in the x-y plane, with the vertices in each cell ordered clockwise, and we defined r i, j = r i − r j and r i, j = r i, j /|r i, j |. To compute the energy gradients, it is convenient to write this expression using any vertex to make the triangularization where cyclic vertex numbering is used (i.e. N + 1 ≡ 1 and −1 ≡ N). Then, Also, the perimeter and its gradient with respect to the position of the vertex i of the same polygon are given by Finally, the different terms of Eq. (B2) are where Eqs. (B8), (B9), and (B10) consider a sum over the three cells at which the vertex i belongs to, and i c + 1 and i c − 1 refer to the next and previous vertex to i, in clockwise counting, belonging to cell c.

Appendix C: Short and long time scales
By performing a simple dimensional analysis we can obtain the relevant time scales of the dynamics, and define useful short time and long time values, τ s and τ l , respectively. The first one allows us to detect the beginning of the instability, while the second allows the non-linear terms, which saturate the eventual instabilities, to act.
We analyze the energy of a single hexagonal cell of equilibrium side a 0 . At time t = 0 it is deformed isotropically such that the new side is a = a 0 + a 1 , with a 1 a 0 . The area (equilibrium area) and perimeter (equilibrium perimeter) are 3 √ 3a 2 /2 (3 √ 3a 2 0 /2) and 6a (6a 0 ), respectively. To simplify, we consider J = 0, in which case the energy of the cell is According to the dynamics of the vertex model, the cell side evolves aṡ where we defined τ A = 1/(K A a 2 0 ) and τ P = 1/(K P ). With the selection of units such that K A = a 0 = 1, we have that τ A = 1 and τ P = 1/p, which is of order 1. Hence, Obviously, for a confluent tissue, the linear and non-linear terms change, and there are parameters for which the coefficients change sign and tissue is stable. Nevertheless, the present analysis allows us to extract the relaxation time scales. The shortest gives the linear evolution, τ 1 ≈ 0.016, and the other two describe the non-linear terms τ 2 ≈ 0.025 and τ 3 ≈ 0.074. If we consider the short time τ s = 0.025, the unstable modes will have grown exponentially, allowing us to identify their effect in the form a change in ellipticity. For the long time τ l = 0.5, the non-linear terms have played a role and the system could have reach a steady state if the non-linear terms saturate the instability.

Appendix D: Comparison between regular and irregular tissues
To compare the dynamics of regular and irregular tissues, we performed simulations for both cases. The results for target perimeter activity, with λ P = +1/2 (predicted line: j = 6p) and λ P = −1/2 (predicted line: j = −6p), can be seen in Figs. 5 and 6, respectively. Although the detailed geometry of the cells change, the flattening parameter and the measure of non-convexity agree remarkable well between regular and irregular tissues, showing that the long wavelength approximation is valid. From Figs. 5b and 6b it is seen that λ P = −1/2 achieves lower values for the standard deviation of the flattening parameter, which results in more rounded cells [ Fig. 6f(II) versus Fig. 5f(II)].  I  II  III  I  II  III FIG. 6. Tissue instabilities obtained in simulations of N = 3000 hexagonal cells with 10% of Gaussian noise over the regular positions (left) and irregular cells (right) (three different tissues considered), under cell target perimeter activity, with λ P = −1/2 and λ A = 0. Same representation as in Fig. 2.

Appendix E: Examples of non-diagonal µ-matrices
Using the expressions in Appendix A it is possible to derive the µ-matrix for different cases. Here, we present some examples where the resulting matrix is non-diagonal, needing the analysis described in Section IV to determine the unstable modes.
The transition line is given by j = 0.569p. Simulation results for irregular tissues can be seen in Fig. 3.
The associated matrix e ps is obtained [Eqs. (10), (11), and (12)] and we compute the curve in parameter space where the minimum eigenvalue of e ps changes its sign. Equivalently we search when the determinant vanishes, finding the linear relation j = −0.583p. Note that, although the Λ and µ matrices are similar to the previous case, the transition line is radically different. Simulation results for irregular tissues can be seen in Fig. 4.
Finally, for a tissue subject to a pure shear pre-deformation, Λ = 1 0.5 0. The line at which the minimum eigenvalue of e dev changes its sign is given by j = −0.769p. Simulation results for irregular tissues can be seen in Fig. 4.