# On the existence of $3$–way $k$–homogeneous Latin trades

Behrooz Bagheri Gh.a,b, Diane Donovanc, and E. S. Mahmoodianb

###### Abstract

A $\mu$-way Latin trade of volume $s$ is a collection of $\mu$ partial Latin
squares $T_{1},T_{2},\ldots,T_{\mu}$, containing exactly the same $s$ filled
cells, such that if cell $(i,j)$ is filled, it contains a different entry in
each of the $\mu$ partial Latin squares, and such that row $i$ in each of the
$\mu$ partial Latin squares contains, set-wise, the same symbols and column
$j$, likewise. It is called $\mu$–way $k$–homogeneous Latin trade, if in each
row and each column $T_{r}$, for $1\leq r\leq\mu,$ contains exactly $k$
elements, and each element appears in $T_{r}$ exactly $k$ times. It is also
denoted by $(\mu,k,m)$ Latin trade, where $m$ is the size of partial Latin
squares.

We introduce some general constructions for $\mu$–way $k$–homogeneous Latin
trades and specifically show that for all $k\leq m$, $6\leq k\leq 13$ and
$k=15$, and for all $k\leq m$, $k=4,\ 5$ (except for four specific values), a
$3$–way $k$–homogeneous Latin trade of volume $km$ exists. We also show that
there are no $(3,4,6)$ Latin trade and $(3,4,7)$ Latin trade. Finally we
present general results on the existence of $3$–way $k$–homogeneous Latin
trades for some modulo classes of $m$.

$a$ Department of Mathematical Sciences

Isfahan University of Technology

84156-83111, Isfahan, I. R. Iran

$b$ Department of Mathematical Sciences

Sharif University of Technology

P. O. Box 11155–9415, Tehran, I. R. Iran

$c$ Department of Mathematics

The University of Queensland

Brisbane 4072

Australia

AMS Subject Classification: 05B15
Keywords: Latin square; Latin trade; $\mu$–way Latin trade; $\mu$–way
$k$–homogeneous Latin trade.

##  1 Introduction

A Latin square $L$ of order $n$ is an $n\times n$ array usually on the set
$N=\\{1,\ldots,n\\}$ where each element of $N$ appears exactly once in each
row and exactly once in each column. We can represent each Latin square as a
subset of $N\times N\times N$,

$L=\\{(i,j;k)\mid\mbox{ element $k$ is located in position }(i,j)\\}.$

A partial Latin square $P$ of order $n$ is an $n\times n$ array of elements
from the set $N$, where each element of $N$ appears at most once in each row
and at most once in each column. The set $S_{P}=\\{(i,j)\mid(i,j;k)\in P\\}$
of the partial Latin square $P$ is called the shape  of $P$ and $|S_{P}|$ is
called the volume of $P$. By ${\cal R}_{P}^{i}$ and ${\cal C}_{P}^{j}$ we mean
the set of entries in row $i$ and column $j$, respectively of $P$. A $\mu$-way
Latin trade, $(T_{1},T_{2},\ldots,T_{\mu})$, of volume $s$ is a collection of
$\mu$ partial Latin squares $T_{1},T_{2},\ldots,T_{\mu}$, containing exactly
the same $s$ filled cells, such that if cell $(i,j)$ is filled, it contains a
different entry in each of the $\mu$ partial Latin squares, and such that row
$i$ in each of the $\mu$ partial Latin squares contains, set-wise, the same
symbols and column $j$, likewise. If $\mu=2$, $(T_{1},T_{2})$ is called a
Latin bitrade. The study of Latin trades and combinatorial trades in general,
has generated much interest in recent years. For a survey on the topic see
[3], [9], and [6].

A $\mu$–way Latin trade which is obtained from another one by deleting its
empty rows and empty columns, is called a $\mu$–way $k$–homogeneous Latin
trade $(\mu\leq k)$ or briefly a $(\mu,k,m)\mbox{ \sf Latin trade}$, if it has
$m$ rows and in each row and each column $T_{r}$, for $1\leq r\leq\mu,$
contains exactly $k$ elements, and each element appears in $T_{r}$ exactly $k$
times.

In Figure 1($a$) a $(3,5,7)\mbox{ Latin trade}$ is demonstrated. The elements
of $T_{2}$ and $T_{3}$ are written as subscripts in the same array as $T_{1}$.
($\bullet$ means the cell is empty.)

${1_{{\bf\displaystyle 2}_{\displaystyle 3}}}$ | ${3_{{\bf\displaystyle 5}_{\displaystyle 2}}}$ | ${5_{{\bf\displaystyle 3}_{\displaystyle 7}}}$ | ${7_{{\bf\displaystyle 1}_{\displaystyle 5}}}$ | ${2_{{\bf\displaystyle 7}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$
---|---|---|---|---|---|---
${{}_{{}_{\displaystyle\bullet}}}$ | ${2_{{\bf\displaystyle 3}_{\displaystyle 4}}}$ | ${4_{{\bf\displaystyle 6}_{\displaystyle 3}}}$ | ${6_{{\bf\displaystyle 4}_{\displaystyle 1}}}$ | ${1_{{\bf\displaystyle 2}_{\displaystyle 6}}}$ | ${3_{{\bf\displaystyle 1}_{\displaystyle 2}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$
${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${3_{{\bf\displaystyle 4}_{\displaystyle 5}}}$ | ${5_{{\bf\displaystyle 7}_{\displaystyle 4}}}$ | ${7_{{\bf\displaystyle 5}_{\displaystyle 2}}}$ | ${2_{{\bf\displaystyle 3}_{\displaystyle 7}}}$ | ${4_{{\bf\displaystyle 2}_{\displaystyle 3}}}$
${5_{{\bf\displaystyle 3}_{\displaystyle 4}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${4_{{\bf\displaystyle 5}_{\displaystyle 6}}}$ | ${6_{{\bf\displaystyle 1}_{\displaystyle 5}}}$ | ${1_{{\bf\displaystyle 6}_{\displaystyle 3}}}$ | ${3_{{\bf\displaystyle 4}_{\displaystyle 1}}}$
${4_{{\bf\displaystyle 5}_{\displaystyle 2}}}$ | ${6_{{\bf\displaystyle 4}_{\displaystyle 5}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${5_{{\bf\displaystyle 6}_{\displaystyle 7}}}$ | ${7_{{\bf\displaystyle 2}_{\displaystyle 6}}}$ | ${2_{{\bf\displaystyle 7}_{\displaystyle 4}}}$
${3_{{\bf\displaystyle 1}_{\displaystyle 5}}}$ | ${5_{{\bf\displaystyle 6}_{\displaystyle 3}}}$ | ${7_{{\bf\displaystyle 5}_{\displaystyle 6}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${6_{{\bf\displaystyle 7}_{\displaystyle 1}}}$ | ${1_{{\bf\displaystyle 3}_{\displaystyle 7}}}$
${2_{{\bf\displaystyle 4}_{\displaystyle 1}}}$ | ${4_{{\bf\displaystyle 2}_{\displaystyle 6}}}$ | ${6_{{\bf\displaystyle 7}_{\displaystyle 4}}}$ | ${1_{{\bf\displaystyle 6}_{\displaystyle 7}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${7_{{\bf\displaystyle 1}_{\displaystyle 2}}}$
${1_{{\bf\displaystyle 2}_{\displaystyle 3}}}$ | ${3_{{\bf\displaystyle 5}_{\displaystyle 2}}}$ | ${5_{{\bf\displaystyle 3}_{\displaystyle 7}}}$ | ${7_{{\bf\displaystyle 1}_{\displaystyle 5}}}$ | ${2_{{\bf\displaystyle 7}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$
---|---|---|---|---|---|---
${{}_{{}_{\displaystyle\bullet}}}$ | $\searrow$ | $\searrow$ | $\searrow$ | $\searrow$ | $\searrow$ | ${{}_{{}_{\displaystyle\bullet}}}$
${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ |  |  |  |  |
| ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ |  |  |  |
|  | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ |  |  |
|  |  | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ |  |
|  |  |  | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ |

($a$) ($b$)

Figure 1: A $(3,5,7)\mbox{ Latin trade}$ and its base row

A $(\mu,k,m)\mbox{ Latin trade}$ $(T_{1},T_{2},\ldots,T_{\mu})$ is called
circulant, if it can be obtained from the elements of its first row, called
the base row and denoted by $\mu$–$B_{m}^{k}$, by permuting the coordinates
cyclically along the diagonals. For example in Figure 1($b$), a
$3$–$B_{7}^{5}$ base row,
$\\{(1,2,3)_{1},(3,5,2)_{2},(5,3,7)_{3},(7,1,5)_{4},(2,7,1)_{5}\\}$, is shown.
Actually if a base row $B=\\{(a_{1},a_{2},\ldots,a_{\mu})_{c_{l}}\mid 1\leq
l\leq k\\}$, where $a_{r}$ and $c_{l}$ $\in\\{1,2,\ldots,m\\},$ is given, we
construct a set of $\mu$ partial Latin squares as in the following manner:

$1\leq r\leq\mu,\ T_{r}=\\{(1+i,c_{l}+i;a_{r}+i)({\rm mod}\ m)|0\leq i\leq
m-1,1\leq l\leq k\\}.$

###### Algorithm 1

To check that $B=\\{(a_{1},a_{2},\ldots,a_{\mu})_{c_{l}}\mid 1\leq l\leq
k\\}$, where $a_{r}$ and $c_{l}$ $\in\\{1,2,\ldots,m\\},$ is a base row of a
$(\mu,k,m)\mbox{ Latin trade}$:
we note that for each $r$, $1\leq r\leq\mu$, ${\cal
R}_{T_{r}}^{1}=\\{a_{r}\mid(a_{1},a_{2},\ldots,a_{\mu})_{c_{l}}\in B\ {\rm
and}\ 1\leq l\leq k\\}$ and ${\cal C}_{T_{r}}^{m}=\\{a_{r}+m-c_{l}\ \equiv
a_{r}-c_{l}({\rm mod}\ m)\mid(a_{1},a_{2},\ldots,a_{\mu})_{c_{l}}\in B\ {\rm
and}\ 1\leq l\leq k\\}$. Now if $B$ satisfies the following conditions, then
it will suffice to be a base row of a $(\mu,k,m)\mbox{ Latin trade}$.

(i) $a_{r}$’s are distinct, for each $(a_{1},a_{2},\ldots,a_{\mu})_{c_{l}}\in
B$.

(ii) $c_{l}$’s are distinct.

(iii) ${\cal R}_{T_{1}}^{1}={\cal R}_{T_{2}}^{1}=\cdots={\cal
R}_{T_{\mu}}^{1}$.

(iv) ${\cal C}_{T_{1}}^{m}={\cal C}_{T_{2}}^{m}=\cdots={\cal
C}_{T_{\mu}}^{m}$.

###### Lemma 1

For each $k\geq\mu$, a $(\mu,k,k)\mbox{ Latin trade}$ exists.

###### Proof.

By taking a Latin square of order $k$ and permuting its rows, cyclically,
$\mu$ times we obtain the desired Latin trade.

A $(\mu,\mu,\mu)$ Latin trade is called a $\mu$–intercalate.

${1_{{\bf\displaystyle 3}_{\displaystyle 2}}}$ | ${2_{{\bf\displaystyle 1}_{\displaystyle 3}}}$ | ${3_{{\bf\displaystyle 2}_{\displaystyle 1}}}$
---|---|---
${3_{{\bf\displaystyle 2}_{\displaystyle 1}}}$ | ${1_{{\bf\displaystyle 3}_{\displaystyle 2}}}$ | ${2_{{\bf\displaystyle 1}_{\displaystyle 3}}}$
${2_{{\bf\displaystyle 1}_{\displaystyle 3}}}$ | ${3_{{\bf\displaystyle 2}_{\displaystyle 1}}}$ | ${1_{{\bf\displaystyle 3}_{\displaystyle 2}}}$

Figure 2: A $3$–intercalate

The following question is of interest.

###### Question 1

For given $m$ and $k$, $m\geq k\geq\mu$, does there exist a $(\mu,k,m)\mbox{
Latin trade}$?

For Latin bitrades, Question 1 is discussed and is answered completely in [4],
[5], [2], [1], and [7]. In this paper applying earlier results we introduce
some general constructions for $(\mu,k,m)\mbox{ Latin trade}$s and
specifically concentrate on the case of $\mu=3$. Our main result is stated in
the following theorem.

###### Theorem 1

All $(3,k,m)\mbox{ Latin trade}$s $(m\geq k\geq 3)$ exist, for

  * •

$k=4$, except for $m=6$ and $7$ and possibly for $m=11$,

  * •

$k=5$, except possibly for $m=6$,

  * •

$6\leq k\leq 13$,

  * •

$k=15$,

  * •

$k\geq 4$ and $m\geq k^{2}$,

  * •

$m$ a multiple of $5$, except possibly for $m=30$,

  * •

$m$ a multiple of $7$, except possibly for $m=42$ and $(3,4,7)\mbox{ Latin
trade}$.

##  2 General constructions

###### Theorem 2

If $l\neq 2,6$ and for each $k\in\\{k_{1},\ldots,k_{l}\\}$ there exists a
$(\mu,k,p)$ Latin trade, then a $(\mu,k_{1}+\cdots+k_{l},lp)\mbox{ Latin
trade}$ exists. (Some $k_{i}$s can possibly be zero.)

###### Proof.

Since $l\neq 2,6$, there exist two $l\times l$ orthogonal Latin squares.
Denote these Latin squares by $L_{1}$ and $L_{2}$, with elements chosen from
the sets $\\{e_{1},e_{2},\ldots,e_{l}\\}$ and
$\\{f_{1},f_{2},\ldots,f_{l}\\}$, respectively. Assume that $L^{*}$ is a
square that is formed by superposing $L_{1}$ and $L_{2}$. We replace each
$(e_{i},f_{j})$ in $L^{*}$ with a $(\mu,k_{j},p)\mbox{ Latin trade}$ whose
elements are from the set $\\{(i-1)p+1,(i-1)p+2,\ldots,ip\\}$. As a result we
obtain a $(\mu,k_{1}+\cdots+k_{l},lp)\mbox{ Latin trade}$.

###### Theorem 3

If the number of mutually orthogonal Latin squares of order $k+1$, ${\rm
MOLS}(k+1)$, is greater than or equal to $\mu+1$, then there exists a
$(\mu,k,k+1)$ Latin trade.

###### Proof.

By Exercise $5.2.11$ of [10] page $103$, there are $\mu$ idempotent ${\rm
MOLS}(k+1).$ If in each of those ${\rm MOLS}$ we delete the main diagonals, we
obtain a $(\mu,k,k+1)$ Latin trade.

Actually by applying results of existence of idempotent ${\rm MOLS}(n)$ ([8],
Section $3.6$, Table $3.83$), we can improve Theorem 3 for the case $\mu=3$ as
follows.

###### Theorem 4

If $k\geq 11$, then there exists a $(3,k,k+1)\mbox{ Latin trade}$.

###### Theorem 5

Any $(\mu,\mu,m)\mbox{ Latin trade}$, ${\bf T}=(T_{1},T_{2},\ldots,T_{\mu})$,
can be partitioned into disjoint $\mu$–intercalates.

###### Proof.

We prove this result by induction. Without loss of generality, let $(1,1;r)\in
T_{r}$ for each $1\leq r\leq\mu$. Therefore $\\{1,2,\ldots,\mu\\}\subset{\cal
R}_{T_{r}}^{i}\cap{\cal C}_{T_{r}}^{i}$ for each $1\leq i,r\leq\mu$. Since
$|{\cal R}_{T_{r}}^{i}|=|{\cal C}_{T_{r}}^{i}|=\mu$ for each $1\leq
i,r\leq\mu$, ${\cal R}_{T_{r}}^{i}={\cal C}_{T_{r}}^{i}=\\{1,2,\ldots,\mu\\}$
for each $1\leq i,r\leq\mu$. Again without loss of generality, let $(i,1;i)\in
T_{1}$ and $(1,j;j)\in T_{1}$ for $1\leq i,j\leq\mu$. This implies that
$\\{(i,j)\mid 1\leq i,j\leq\mu\\}$ is a subset of shape of $T_{1}$. Therefore
subarray $\\{(i,j)\mid 1\leq i,j\leq\mu\\}$ with elements
$\\{1,2,\ldots,\mu\\}$ is a $\mu$–intercalate. We can apply the same argument
to the $(m-\mu)\times(m-\mu)$ subsquare obtained by removing rows
$1,2,\ldots,\mu$ and columns $1,2,\ldots,\mu$. This completes the proof.

###### Corollary 1

For every $m\geq 1$, there exists a $(\mu,k,m)\mbox{ Latin trade}$ with
$k=\mu$, if and only if $k|m$.

###### Theorem 6

Assume that $m_{i}\geq k_{i}$, for $i=1,2$. If there exists a
$(\mu_{i},k_{i},m_{i})$ Latin trade for $i=1,2$, then there exists a
$(\mu_{1}\mu_{2},k_{1}k_{2},m_{1}m_{2})\mbox{ Latin trade}$.

###### Proof.

We construct a $(\mu_{1}\mu_{2},k_{1}k_{2},m_{1}m_{2})\mbox{ Latin trade}$ in
the following way:
Suppose $(T_{1},T_{2},\ldots,T_{\mu_{1}})$ is a $(\mu_{1},k_{1},m_{1})\mbox{
Latin trade}$ and ${\bf U}=(U_{1},U_{2},\ldots,U_{\mu_{2}})$ is a
$(\mu_{2},k_{2},m_{2})\mbox{ Latin trade}$. For each entry $i$ in
$T_{1},T_{2},\ldots,T_{\mu_{1}}$, we replace $i$ with a copy of ${\bf U}$
where elements are chosen from the set
$\\{(i-1)m_{2}+1,(i-1)m_{2}+2,\ldots,im_{2}\\}$; replace the empty cells in
$T_{1},T_{2},\ldots,T_{\mu_{1}}$ with an empty $m_{2}\times m_{2}$ array. As a
result we obtain a $(\mu_{1}\mu_{2},k_{1}k_{2},m_{1}m_{2})\mbox{ Latin
trade}$.

###### Corollary 2

Suppose $k=k_{1}k_{2}$ and $m=m_{1}m_{2}$ where $m_{i}\geq k_{i}\geq 2$, for
$i=1,2$. Then there exists a $(4,k,m)\mbox{ Latin trade}$, provided that if
$k_{j}=2$, for some $j$, then $m_{j}$ must be assumed to be even.

###### Proof.

It is shown that Latin homogeneous bitrades (i.e $(2,k,m)\mbox{ Latin trade}$)
exist for all $m\geq k\geq 3$ and for all even $m$, when $k=2$. (See [4], [5],
[2], [1], and [7].)

###### Theorem 7

For every $k$, if there exists a $(\mu,k,m)\mbox{ Latin trade}$ and a
$(\mu,k,n)$ Latin trade, then there exists a $(\mu,k,m+n)\mbox{ Latin trade}$.

###### Proof.

Let ${\bf T_{1}}$ be a $(\mu,k,m)\mbox{ Latin trade}$ and ${\bf T_{2}}$ be a
$(\mu,k,n)\mbox{ Latin trade}$ such that the elements of ${\bf T_{1}}$ are in
the set $\\{1,\ldots,m\\}$ and the elements of ${\bf T_{2}}$ are chosen from
the set $\\{m+1,\ldots,m+n\\}$. Therefore, the following Latin trade is a
$(\mu,k,m+n)\mbox{ Latin trade}$.

$\begin{array}[]{|@{\hspace{2pt}}c@{\hspace{2pt}}@{\hspace{1pt}}c@{\hspace{1pt}}
@{\hspace{1pt}}c@{\hspace{1pt}}@{\hspace{1pt}}c@{\hspace{1pt}}
|}\hline\cr\hskip 2.0pt\lx@intercol\hfil{\bf T_{1}}\hfil\hskip 2.0pt\hskip
1.0&\hfil\hskip 1.0pt\hskip 1.0&\hfil\hskip 1.0pt\hskip 1.0&\hfil\hskip 1.0\\\
\hskip 2.0pt\lx@intercol\hfil\hfil\hskip 2.0pt\hskip 1.0&{\bf
T_{2}}\hfil\hskip 1.0pt\hskip 1.0&\hfil\hskip 1.0pt\hskip 1.0&\hfil\hskip
1.0\\\ \hline\cr\end{array}$

###### Corollary 3

If the number of ${\rm MOLS}(k+1)\geq\mu+1$, then for each $m$ where $m\geq
k^{2}$, there exists a $(\mu,k,m)\mbox{ Latin trade}$.

###### Proof.

If $m\geq k^{2}$, then we can write $m$ as $m=rk+s(k+1)$, where $r,s\geq 0$.
Theorem 7 and Theorem 3 lead us to a conclusion.

By Theorems 7 and 4 we have:

###### Corollary 4

If $k\geq 11$, then for each $m$ where $m\geq k^{2}$, there exists a $(3,k,m)$
Latin trade.

###### Theorem 8

Consider an arbitrary natural number $k$. If for every ${k+1}\leq l\leq 2k-1$
there exists a $(\mu,k,l)\mbox{ Latin trade}$, then for any $m\geq k$ there
exists a $(\mu,k,m)$ Latin trade.

###### Proof.

For every $m\geq 2k$, we can write $m=rk+sl$, where $r,s\geq 0$ and ${k+1}\leq
l\leq 2k-1$. Since there exist a $(\mu,k,k)\mbox{ Latin trade}$ and a
$(\mu,k,l)\mbox{ Latin trade}$, by Theorem 7 we conclude that there exists a
$(\mu,k,m)$ Latin trade.

##  3 $\mu=3$

In this section we apply the above constructions to establish the existence of
$3$–way $k$–homogeneous Latin trades for specific values of $k$, and when $m$
is a multiple of $5$ or $7$. We also show that there is no $(3,4,6)$ Latin
trade.

###  3.1 Small even $k$

###### Proposition 1

There exists a $(3,4,m)\mbox{ Latin trade}$ for every $m\geq 4$, except
possibly for $m=6,7$ and $11$.

###### Proof.

By Lemma 1 and Theorem 3 there exist a $(3,4,4)\mbox{ Latin trade}$ and a
$(3,4,5)$ Latin trade, respectively. Since $8=2\times 4,\ 9=4+5,\ 10=2\times
5,\ 12=3\times 4,\ 13=2\times 4+5,\ 14=4+2\times 5,$ and $15=3\times 5$;
Theorem 7 results that there exist $(3,4,m)\mbox{ Latin trade}$s for
$m=8,9,10,12,13,14$, and $15$. Since the number ${\rm MOLS}(5)=4$, then by
Corollary 3 there exists a $(3,4,m)$ Latin trade, for every $m\geq 16$.

###### Proposition 2

There is no $(3,4,6)$ Latin trade.

###### Proof.

By contradiction. Suppose $T=(T_{1},T_{2},T_{3})$ is a $(3,4,6)\mbox{ Latin
trade}$. By applying some permutations on rows and columns, if necessary, we
may assume that all cells containing the element 1 form a $4\times 4$ array
minus a transversal $\tau$, which will be labeled $L$. For example in Figure 3
one of the possible positions of 1 is shown. Note that there are 12 cells in
$L$ each of which has a 1 in one of the $T_{i}$’s. In what follows the
argument is based only on the assumption that in each of those cells there
exists one 1 from one of the $T_{i}$’s. ($\bullet$ means the cell is empty.)

${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ |  |
---|---|---|---|---|---
${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ |  |
${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ |  |
${{}_{{}_{\displaystyle\bullet}}}$ | ${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ |  |
|  |  |  |  |
|  |  |  |  |

Figure 3: Positions of 1 in $T=(T_{1},T_{2},T_{3})$

In the first stage we show that the cells of $\tau$ in $T$ are empty. Suppose
without loss of generality the cell $T_{14}$ in $\tau$ is not empty. Then
$T_{54}$ and $T_{64}$ must be empty. Thus at least 4 cells of
$\\{T_{51},T_{52},T_{53},T_{61},T_{62},T_{63}\\}$ must be filled. Then by
pigeonhole principal there exists a column in $T$ with at least 5 filled
cells, a contradiction. So all cells of $\tau$ are empty. Therefore exactly 4
cells of $\\{T_{51},T_{52},T_{53},T_{54},T_{61},T_{62},T_{63},T_{64}\\}$ are
filled, and from $T$ being $4$-homogeneous all the cells:
$\\{T_{55},T_{56},T_{65},T_{66}\\}$ are filled.

In the second stage we show that no element, other than 1, appears more than
two times in any row or in any column of $L$. For example let us denote by
$\\{1,x,y,z\\}$, the elements which appear in the first row and without loss
of generality $T_{15}$ is another filled cell of that row. In contrary, assume
that $x$ appears three times in the first row of $L$, i.e. in the cells
$T_{11},T_{12},$ and $T_{13}$. This leaves only two elements $y$ and $z$ to
appear in $T_{15}$, which is a contradiction for $T$ being a $3$-way Latin
trade. So each of the elements other than 1, either does not appear in a row
of $L$ or it appears exactly two times in a row of $L$. Now each element other
than 1 if it appears in $L$, it occupies 4, 6, or 8 cells.

In the third stage we show that no element occupies 6 or 8 cells of $L$. If an
element, say $u\neq 1$ appears 8 times in $L$, then since $u$ appears 2 times
in each row and in each column of $L$, so it appears once in each row of the
$[1,\ldots,4]\times[5,6]$ block. This means that $u$ appears at least 16 times
in $T$, which is a contradiction. If $u\neq 1$ appears 6 times in $L$ then
three rows and three columns of $L$ each contains $u$ twice. So without loss
of generality one of the following cases happens.

${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | u | ${{}_{{}_{\displaystyle\bullet}}}$
---|---|---|---|---|---
${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | u | ${{}_{{}_{\displaystyle\bullet}}}$
${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | u | ${{}_{{}_{\displaystyle\bullet}}}$
${{}_{{}_{\displaystyle\bullet}}}$ | ${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ |  |
|  |  |  |  |
|  |  |  |  |
${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | u | ${{}_{{}_{\displaystyle\bullet}}}$
---|---|---|---|---|---
${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | u | ${{}_{{}_{\displaystyle\bullet}}}$
${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | u
${{}_{{}_{\displaystyle\bullet}}}$ | ${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ |  |
|  |  |  |  |
|  |  |  |  |

($a$) ($b$)

Figure 4: Positions of u in the fifth and sixth columns of $T$

In case ($a$) the fifth column has at least 5 filled cells which is a
contradiction. In case ($b$) there are five columns of $T$ which have $u$ and
since each column containing an $u$ will contain 3 of them, so there are at
least 15 cells containing $u$ in $T$, which is a contradiction.

Now we have shown that each $u\neq 1$ if it appears in $L$, it appears exactly
4 times. The array $L$ has exactly $36-12=24$ places for elements different
from 1 to occupy while the 5 other elements can fill at most $5\times 4=20$
places, which is a contradiction.

###### Proposition 3

There is no $(3,4,7)$ Latin trade.

###### Proof.

By contradiction. Suppose $T=(T_{1},T_{2},T_{3})$ is a $(3,4,7)\mbox{ Latin
trade}$. By applying some permutations on rows and columns, if necessary, we
may assume that all cells containing the element 1 form a $4\times 4$ array
minus a transversal $\tau$, which will be labeled $L$. For example in Figure 5
one of the possible positions of 1 is shown. Note that there are 12 cells in
$L$ each of which has a 1 in one of the $T_{i}$’s. In what follows the
argument is based only on the assumption that in each of those cells there
exists one 1 from one of the $T_{i}$’s. ($\bullet$ means the cell is empty.)

${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ |  |  |
---|---|---|---|---|---|---
${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ |  |  |
${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ |  |  |
${{}_{{}_{\displaystyle\bullet}}}$ | ${1_{{\bf\displaystyle.}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle 1}_{\displaystyle.}}}$ | ${._{{\bf\displaystyle.}_{\displaystyle 1}}}$ |  |  |
|  |  |  |  |  |
|  |  |  |  |  |
|  |  |  |  |  |

Figure 5: Positions of 1 in $T=(T_{1},T_{2},T_{3})$

If we focus on the placement of the remaining filled cells in $T$, we see that
rows $1$ to $4$ of $T$ each have one additional filled cell in one of columns
$5,6$ or $7$. Likewise for columns $1$ to $4$ of rows $5,6$ or $7$. Further,
the subsquare defined by the intersection of rows $5,6,$ and $7$ with columns
$5,6,$ and $7$, can have at most three filled cells in any row or column.
Hence it follows that without loss of generality columns $5$ has two filled
cell in rows $1$ to $4$ (similarly row $5$ has two filled cells in columns $1$
to $4$) and columns $6$ and $7$ have one filled cell in rows $1$ to $4$
(similarly rows $6$ and $7$ have one filled cell in columns $1$ to $4$). Thus
we may assume cell $(5,5)$ is empty and one possible distribution of empty
cells (one out of $36$) is:

|  |  | ${{}_{{}_{\displaystyle\bullet}}}$ |  | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$
---|---|---|---|---|---|---
|  | ${{}_{{}_{\displaystyle\bullet}}}$ |  |  | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$
| ${{}_{{}_{\displaystyle\bullet}}}$ |  |  | ${{}_{{}_{\displaystyle\bullet}}}$ |  | ${{}_{{}_{\displaystyle\bullet}}}$
${{}_{{}_{\displaystyle\bullet}}}$ |  |  |  | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ |
|  | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ |  |
${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ |  |  |  |
${{}_{{}_{\displaystyle\bullet}}}$ | ${{}_{{}_{\displaystyle\bullet}}}$ |  | ${{}_{{}_{\displaystyle\bullet}}}$ |  |  |

We can assume that the cell $T_{15}$ contains symbols $2,3,4$. Then the first
row must contain only symbols $1,2,3,4,$ and these are distributed among the
four filled cells (in the first row) according to one of three possible ways:

123, 124, 134, 234 (or 123, 134, 124, 234)

124, 134, 123, 234 (or 124, 123, 134, 234)

134, 124, 123, 234 (or 134, 123, 124, 234)

The idea is to label the filled columns with one of these configurations, to
label the first row 1234, and then attempt to complete the labeling of the
rows and columns as follows:

$\bullet$ each row and column is labeled by 4 elements from
$\\{1,\ldots,7\\}$,

$\bullet$ the first 4 rows and first 4 columns contain 1 in its label,

$\bullet$ first row is labeled $\\{1,2,3,4\\}$,

$\bullet$ columns with filled cells in the first row are filled as above,

$\bullet$ for any $i$, the number $i$ appears in precisely 4 row labels and in
precisely 4 column labels,

$\bullet$ if the cell $T_{ij}$ is filled, $A$ is the label of row $i$ and $B$
is the label of row $j$, then $|A\cup B|\leq 5$ (because the cell $T_{ij}$
contains three elements of $A\cap B$).

By applying a depth-first search, we found no solutions (Indeed, we tried all
36 distributions of filled cells and all three configurations in the first
row). The search takes a minute with no optimization. So, it is already
impossible to distribute elements in rows and columns according to the
restrictions of the $(3,4,7)\mbox{ Latin trade}$ disregarding how the cell
symbols are distributed among the three components of the purported Latin
trade. Therefore, there is no $(3,4,7)\mbox{ Latin trade}$.

At this point we will show the existence of some $(3,k,m)\mbox{ Latin
trade}$s. For this purpose we will need some small cases. We have found base
rows of those Latin trades computationally, sometimes by trial and errors. But
we have checked all of them by Algorithm 1.

###### Theorem 9

If $k=6,8,10$ and $12$ then there exists a $(3,k,m)\mbox{ Latin trade}$ for
every $m\geq k$.

###### Proof.

We will show for the given $k$, there exist $(\mu,k,l)\mbox{ Latin trade}$s
for $l$, where ${k+1}\leq l\leq 2k-1$. Then by Theorem 8, we will get all
$m\geq k$ where $k=6,8,10$, and $12$.

  * •

$k=6$.

If $8\leq m=2l\leq 10$, by Corollary 2 a $(3,6,m)\mbox{ Latin trade}$ exists.

And the following are the base rows of a $(3,6,m)\mbox{ Latin trade}$ for
$m=7,9,11$:

$3$–$B_{7}^{6}=\\{(1,5,4)_{1},(3,4,2)_{2},(5,3,1)_{3},(7,2,5)_{4},(2,1,7)_{5},(4,7,3)_{6}\\}$,
$3$–$B_{9}^{6}=\\{(1,8,3)_{1},(3,2,1)_{2},(2,5,6)_{3},(6,3,2)_{4},(8,6,5)_{5},(5,1,8)_{7}\\}$,
$3$–$B_{11}^{6}=\\{(1,6,3)_{1},(3,2,7)_{2},(6,4,1)_{3},(2,7,4)_{4},(7,3,6)_{5},(4,1,2)_{10}\\}.$

  * •

$k=8$.

If $10\leq m=2l\leq 14$, by Corollary 2 a $(3,8,m)\mbox{ Latin trade}$ exists.

And the following are the base rows of a $(3,8,m)\mbox{ Latin trade}$ for
$m=9,11,13,15$:

$3$–$B_{9}^{8}=\\{(1,8,7)_{1},(3,2,9)_{2},(2,4,3)_{3},(7,1,6)_{4},(9,7,4)_{5},(8,9,1)_{6},(4,6,8)_{7},$
$(6,3,2)_{8}\\}$,

$3$–$B_{11}^{8}=\\{(1,5,4)_{1},(3,2,11)_{2},(2,4,5)_{3},(6,1,3)_{4},(8,3,2)_{5},(4,8,6)_{6},(11,6,8)_{7},$
$(5,11,1)_{8}\\}$,

$3$–$B_{13}^{8}=\\{(1,5,3)_{1},(3,1,5)_{2},(2,6,11)_{3},(6,4,2)_{4},(8,3,4)_{5},(4,8,6)_{6},(11,2,8)_{7},$
$(5,11,1)_{10}\\}$,

$3$–$B_{15}^{8}=\\{(1,11,4)_{1},(3,2,6)_{2},(2,4,3)_{3},(6,7,2)_{4},(8,3,7)_{5},(4,8,1)_{6},(11,6,8)_{7},$
$(7,1,11)_{12}\\}.$

  * •

$k=10$.

If $12\leq m=2l\leq 18$, by Corollary 2 a $(3,10,m)\mbox{ Latin trade}$
exists.

And the following are the base rows of a $(3,10,m)\mbox{ Latin trade}$ for
$m=13,15,17,19$:

$3$–$B_{13}^{10}=\\{(1,11,6)_{1},(3,2,13)_{2},(2,4,3)_{3},(6,8,7)_{4},(8,7,4)_{5},(4,5,2)_{6},(11,3,8)_{7},$
$(13,6,5)_{8},(5,1,11)_{9},(7,13,1)_{10}\\}$,

$3$–$B_{15}^{10}=\\{(1,6,5)_{1},(3,2,4)_{2},(2,4,14)_{3},(6,8,3)_{4},(8,1,2)_{5},(4,3,6)_{6},(11,5,8)_{7},$
$(5,7,11)_{8},(14,11,7)_{9},(7,14,1)_{11}\\}$,

$3$–$B_{17}^{10}=\\{(1,6,4)_{1},(3,2,6)_{2},(2,7,14)_{3},(6,1,2)_{4},(8,4,5)_{5},(4,8,3)_{6},(11,5,8)_{7},$
$(5,11,7)_{8},(14,3,11)_{9},(7,14,1)_{13}\\}$,

$3$–$B_{19}^{10}=\\{(1,6,2)_{1},(3,2,6)_{2},(2,4,14)_{3},(6,8,7)_{4},(8,7,3)_{5},(4,3,5)_{6},(11,5,4)_{7},$
$(5,11,8)_{8},(14,1,11)_{9},(7,14,1)_{15}\\}.$

  * •

$k=12$.

If $14\leq m=2l\leq 22$ or $m=15,21$, by Corollary 2 a $(3,12,m)\mbox{ Latin
trade}$ exists.

And the following are the base rows of a $(3,12,m)\mbox{ Latin trade}$ for
$m=17,19,23$:

$3$–$B_{17}^{12}=\\{(1,16,4)_{1},(3,7,2)_{2},(2,4,14)_{3},(6,8,3)_{4},(8,5,11)_{5},(4,3,10)_{6},(11,1,8)_{7}$,
$(5,14,6)_{8},(14,11,5)_{9},(16,6,7)_{10},(7,10,16)_{11},(10,2,1)_{16}\\}$,

$3$–$B_{19}^{12}=\\{(1,16,7)_{1},(3,2,6)_{2},(2,4,3)_{3},(6,9,1)_{4},(8,7,4)_{5},(4,3,11)_{6},(11,5,2)_{7},$
$(5,11,9)_{8},(14,8,5)_{9},(16,14,8)_{10},(7,6,14)_{11},(9,1,16)_{14}\\}$,

$3$–$B_{23}^{12}=\\{(1,7,5)_{1},(3,2,8)_{2},(2,4,1)_{3},(6,9,3)_{4},(8,1,7)_{5},(4,3,9)_{6},(11,5,2)_{7},$
$(5,11,4)_{8},(14,8,6)_{9},(16,14,11)_{10},(7,6,16)_{11},(9,16,14)_{14}\\}.$

###  3.2 Small odd $k$

###### Proposition 4

There exists a $(3,5,m)\mbox{ Latin trade}$ for every $m\geq 5$, except
possibly $m=6$.

###### Proof.

By Lemma 1 there exists a $(3,5,5)\mbox{ Latin trade}$. The following are the
base rows of a $(3,5,m)\mbox{ Latin trade}$ for $m=7,8,9,11$:

$3$–$B_{7}^{5}=\\{(1,3,2)_{1},(3,2,5)_{2},(5,7,3)_{3},(7,5,1)_{4},(2,1,7)_{5}\\}$,
$3$–$B_{8}^{5}=\\{(1,6,2)_{1},(3,2,4)_{2},(2,4,3)_{3},(6,3,1)_{4},(4,1,6)_{7}\\}$,
$3$–$B_{9}^{5}=\\{(1,4,3)_{1},(4,3,8)_{2},(7,1,4)_{4},(3,8,7)_{6},(8,7,1)_{7}\\}$,
$3$–$B_{11}^{5}=\\{(1,6,9)_{1},(9,2,11)_{5},(11,1,6)_{6},(2,11,1)_{7},(6,9,2)_{9}\\}.$

By Theorem 7, a $(3,5,10)\mbox{ Latin trade}$ exists. So a $(3,5,m)\mbox{
Latin trade}$ exists for 5 consecutive values $m\in\\{7,8,\ldots,11\\}$. Thus
a $(3,5,m)\mbox{ Latin trade}$ exists for all $m\geq 7$ by Theorem 7.

###### Theorem 10

If $k=7,9,11$ and $13$ then there exists a $(3,k,m)\mbox{ Latin trade}$ for
every $m\geq k$.

###### Proof.

We introduce the following base rows:

  * •

$k=7$.

    * $m\geq 8$:

$3$–$B_{m}^{7}=\\{(1,4,2)_{1},(3,1,4)_{2},(2,3,6)_{3},(6,5,1)_{4},(8,2,3)_{5},(4,8,5)_{6},(5,6,8)_{8}\\}.$

  * •

$k=9$.

    * $m=10$:

$3$–$B_{10}^{9}=\\{(1,7,9)_{1},(3,2,8)_{2},(2,4,5)_{3},(7,6,4)_{4},(9,3,2)_{5},(8,9,7)_{6},(4,1,6)_{7},$
$(6,5,1)_{8},(5,8,3)_{9}\\}.$

    * $m\geq 11$:

$3$–$B_{m}^{9}=\\{(1,5,4)_{1},(3,4,6)_{2},(2,3,1)_{3},(6,2,5)_{4},(8,1,2)_{5},(4,7,8)_{6},(11,6,3)_{7},$
$(5,11,7)_{8},(7,8,11)_{11}\\}.$

  * •

$k=11$.

    * $m\geq 11$:

$3$–$B_{m}^{11}=\\{(6,1,2)_{1},(1,7,4)_{2},(7,2,1)_{3},(2,8,7)_{4},(8,3,10)_{5},(3,9,5)_{6},(9,4,11)_{7},$
$(4,10,3)_{8},(10,5,9)_{9},(5,11,6)_{10},(11,6,8)_{11}\\}.$

  * •

$k=13$.

    * $m\geq 13$:

$3$–$B_{m}^{13}=\\{(7,1,2)_{1},(1,8,4)_{2},(8,2,1)_{3},(2,9,3)_{4},(9,3,8)_{5},(3,10,11)_{6},(10,4,13)_{7},$
$(4,11,12)_{8},(11,5,6)_{9},(5,12,10)_{10},(12,6,5)_{11},(6,13,7)_{12},(13,7,9)_{13}\\}.$

###### Theorem 11

If $k=15$ and $m\geq 15$ then there exists a $(3,15,m)\mbox{ Latin trade}$.

###### Proof.

By Lemma 1, Theorem 3 and Corollary 2, we have a $(3,15,m)\mbox{ Latin trade}$
for $m=15,16,18$ and $20$. The following is a base row of a $(3,15,m)\mbox{
Latin trade}$ for $m\geq 21$:

$3$–$B_{m}^{15}=\\{(1,5,4)_{1},(3,1,2)_{2},(2,9,11)_{3},(6,11,3)_{4},(8,6,7)_{5},(4,14,10)_{6},(11,4,8)_{7},$
$(5,3,6)_{8},(14,7,5)_{9},(16,10,1)_{10},(7,2,16)_{11},(19,8,9)_{12},(21,16,19)_{13},$
$(9,19,21)_{14},(10,21,14)_{19}\\}.$

The following are the base rows of a $(3,15,m)\mbox{ Latin trade}$ for
$m=17,19$:

$3$–$B_{17}^{15}=\\{(5,2,12)_{3},(7,15,11)_{4},(9,17,4)_{5},(11,13,14)_{6},(13,16,5)_{7},(15,11,13)_{8},$
$(17,14,6)_{9},(2,12,16)_{10},(4,9,7)_{11},(6,8,15)_{12},(8,10,17)_{13},$
$(10,7,8)_{14},(12,6,10)_{15},(14,5,9)_{16},(16,4,2)_{17}\\},$

$3$–$B_{19}^{15}=\\{(1,2,11)_{1},(3,4,2)_{2},(5,17,4)_{3},(7,10,9)_{4},(9,15,14)_{5},(11,9,13)_{6},$
$(13,19,10)_{7},(15,13,5)_{8},(17,6,1)_{9},(19,14,3)_{10},(2,11,15)_{11},$
$(4,1,7)_{12},(6,3,19)_{13},(10,7,17)_{15},(14,5,6)_{17}\\}.$

###  3.3 General cases

###### Theorem 12

Let $m\equiv 1({\rm mod\ 6})$ and $m\geq 7$. Then there exists a $(3,m-2,m)$
Latin trade.

###### Proof.

The following is a base row of a $(2,m-2,m)\mbox{ Latin trade}$:

$2$–$B_{m}^{m-2}=\bigcup_{i=0}^{(m-13)/6}\\{(6i+2,6i+3)_{3i+1},(6i+4,6i+2)_{3i+2},(6i+3,6i+4)_{3i+3},$
$(6i+5,6i+6)_{(m+3)/2+3i+1},(6i+7,6i+5)_{(m+3)/2+3i+2},(6i+6,6i+7)_{(m+3)/2+3i+3}\\}$
$\bigcup\\{(m-5,m-4)_{(m-7)/2+1},(m-2,m-5)_{(m-7)/2+2},(m-4,m)_{(m-7)/2+3},$
$(1,m-2)_{(m-7)/2+4},(m,1)_{(m-7)/2+5}\\}.$

Now, for $1\leq i\leq m-2$ we put $2i-1({\rm mod\ m})$ in $i$-th cell of
$2$–$B_{m}^{m-2}$, as a result we obtain a base row of a $(3,m-2,m)\mbox{
Latin trade}$.

###### Example 1

As an example of the previous theorem, the following is a base row of a
$(3,11,13)\mbox{ Latin trade}$:

$3$–$B_{13}^{11}=\\{(1,2,3)_{1},(3,4,2)_{2},(5,3,4)_{3},(7,8,9)_{4},(9,11,8)_{5},(11,9,13)_{6},(13,1,11)_{7},$
$(2,13,1)_{8},(4,5,6)_{9},(6,7,5)_{10},(8,6,7)_{11}\\}.$

###### Theorem 13

For every $m=5l$ and $4\leq k\leq m$, $l\neq 6$, there exists a $(3,k,m)$
Latin trade.

###### Proof.

The theorem trivially holds for $l=1$. If $l=2$, then by Theorem 3, Theorem 7,
Theorem 10 and Theorem 9, we can construct a $(3,k,10)\mbox{ Latin trade}$ for
every $4\leq k\leq 10$. By Theorems 9 and 10 there exists a $(3,k,m)\mbox{
Latin trade}$ for $k=6,7$ and $11$, so suppose that $k\neq 6,7$ and $11$.

We may also assume that $m>k$.
We have the following cases to consider, each case follows from Theorem 2:

  * •

$k=5l^{{}^{\prime}}$.
We set $k_{i}=5$ for $1\leq i\leq l^{{}^{\prime}}$ and $k_{i}=0$ for
$l^{{}^{\prime}}+1\leq i\leq l$ and $p=5$.

  * •

$k=5l^{{}^{\prime}}+1$.
We set $k_{i}=5$ for $1\leq i\leq l^{{}^{\prime}}-3$ and $k_{i}=4$ for
$l^{{}^{\prime}}-2\leq i\leq l^{{}^{\prime}}+1$ and $k_{i}=0$ for
$l^{{}^{\prime}}+2\leq i\leq l$ and $p=5$.

  * •

$k=5l^{{}^{\prime}}+2$.
We set $k_{i}=5$ for $1\leq i\leq l^{{}^{\prime}}-2$ and $k_{i}=4$ for
$l^{{}^{\prime}}-1\leq i\leq l^{{}^{\prime}}+1$ and $k_{i}=0$ for
$l^{{}^{\prime}}+2\leq i\leq l$ and $p=5$.

  * •

$k=5l^{{}^{\prime}}+3$.
We set $k_{i}=5$ for $1\leq i\leq l^{{}^{\prime}}-1$,
$k_{l^{{}^{\prime}}}=k_{l^{{}^{\prime}}+1}=4$, and $k_{i}=0$ for
$l^{{}^{\prime}}+2\leq i\leq l$, and $p=5$.

  * •

$k=5l^{{}^{\prime}}+4$.
We set $k_{i}=5$ for $1\leq i\leq l^{{}^{\prime}}$, $k_{l^{{}^{\prime}}+1}=4$,
and $k_{i}=0$ for $l^{{}^{\prime}}+2\leq i\leq l$, and $p=5$.

###### Theorem 14

For every $m=7l$ and $5\leq k\leq m$, $l\neq 6$, there exists a $(3,k,m)$
Latin trade.

###### Proof.

The theorem trivially holds for $l=1$. If $l=2$, then by Theorem 7, Theorem 10
and Theorem 9, we can construct a $(3,k,14)\mbox{ Latin trade}$ for every
$5\leq k\leq 14$. For $l\neq 2,6$ by Theorems 9 and 10 there exists a
$(3,k,m)\mbox{ Latin trade}$ for $k=8,9$, so suppose that $k\neq 8,9$.

We may also assume that $m>k$.
We have the following cases to consider, each case follows from Theorem 2:

  * •

$k=7l^{{}^{\prime}}$.
We set $k_{i}=7$ for $1\leq i\leq l^{{}^{\prime}}$ and $k_{i}=0$ for
$l^{{}^{\prime}}+1\leq i\leq l$ and $p=7$.

  * •

$k=7l^{{}^{\prime}}+1$.
We set $k_{i}=7$ for $1\leq i\leq l^{{}^{\prime}}-2$ and $k_{i}=5$ for
$l^{{}^{\prime}}-1\leq i\leq l^{{}^{\prime}}+1$ and $k_{i}=0$ for
$l^{{}^{\prime}}+2\leq i\leq l$ and $p=7$.

  * •

$k=7l^{{}^{\prime}}+2$.
We set $k_{i}=7$ for $1\leq i\leq l^{{}^{\prime}}-2$ and
$k_{l^{{}^{\prime}}-1}=k_{l^{{}^{\prime}}}=5$, $k_{l^{{}^{\prime}}+1}=6$ and
$k_{i}=0$ for $l^{{}^{\prime}}+2\leq i\leq l$ and $p=7$.

  * •

$k=7l^{{}^{\prime}}+3$.
We set $k_{i}=7$ for $1\leq i\leq l^{{}^{\prime}}-1$,
$k_{l^{{}^{\prime}}}=k_{l^{{}^{\prime}}+1}=5$, and $k_{i}=0$ for
$l^{{}^{\prime}}+2\leq i\leq l$, and $p=7$.

  * •

$k=7l^{{}^{\prime}}+4$.
We set $k_{i}=7$ for $1\leq i\leq l^{{}^{\prime}}-1$, $k_{l^{{}^{\prime}}}=6$,
$k_{l^{{}^{\prime}}+1}=5$ and $k_{i}=0$ for $l^{{}^{\prime}}+2\leq i\leq l$,
and $p=7$.

  * •

$k=7l^{{}^{\prime}}+5$.
We set $k_{i}=7$ for $1\leq i\leq l^{{}^{\prime}}$, $k_{l^{{}^{\prime}}+1}=5$,
and $k_{i}=0$ for $l^{{}^{\prime}}+2\leq i\leq l$, and $p=7$.

  * •

$k=7l^{{}^{\prime}}+6$.
We set $k_{i}=7$ for $1\leq i\leq l^{{}^{\prime}}$, $k_{l^{{}^{\prime}}+1}=6$,
and $k_{i}=0$ for $l^{{}^{\prime}}+2\leq i\leq l$, and $p=7.$

Now by the results given above we have proved Theorem 1, given at the end of
the Introduction.

Acknowledgement. We appreciate very useful comments of anonymous referee, who
also added the proof of Proposition 3. Also we thank to Amir Hooshang
Hosseinpoor and Mayssam Mohammadi Nevisi for their computer programming.

## References

  * [1] Behrooz Bagheri Gh. and E. S. Mahmoodian.  On the existence of $k$–homogeneous Latin bitrades.  Util. Math., 85:333–345, 2011.
  * [2] Richard Bean, Hoda Bidkhori, Maryam Khosravi, and E. S. Mahmoodian.  $k$-homogeneous Latin trades.  Bayreuth. Math. Schr., 74:7–18, 2005.
  * [3] Elizabeth J. Billington.  Combinatorial trades: a survey of recent results.  In Designs, 2002, volume 563 of Math. Appl., pages 47–67. Kluwer Acad. Publ., Boston, MA, 2003.
  * [4] Nicholas Cavenagh, Diane Donovan, and Aleš Drápal.  3-homogeneous Latin trades.  Discrete Math., 300(1-3):57–70, 2005.
  * [5] Nicholas Cavenagh, Diane Donovan, and Aleš Drápal.  4-homogeneous Latin trades.  Australas. J. Combin., 32:285–303, 2005.
  * [6] Nicholas J. Cavenagh.  The theory and application of Latin bitrades: a survey.  Math. Slovaca, 58(6):691–718, 2008.
  * [7] Nicholas J. Cavenagh and Ian M. Wanless.  On the number of transversals in Cayley tables of cyclic groups.  Discrete Appl. Math., 158(2):136–146, 2010.
  * [8] Charles J. Colbourn and Jeffrey H. Dinitz, editors.  Handbook of combinatorial designs.  Discrete Mathematics and its Applications (Boca Raton). Chapman & Hall/CRC, Boca Raton, FL, second edition, 2007.
  * [9] A. D. Keedwell.  Critical sets in Latin squares and related matters: an update.  Util. Math., 65:97–131, 2004.
  * [10] C. C. Lindner and C. A. Rodger.  Design theory.  CRC Press LLC, 1997.


# Modelling realistic 3D deformations of simple epithelia in dynamic
homeostasis

Domenic P.J. Germano  Stuart T. Johnston  Edmund J. Crampin  James M. Osborne
<EMAIL_ADDRESS>School of Mathematics and Statistics, University of
Melbourne, Parkville, Victoria 3010, Australia  Systems Biology Laboratory,
School of Mathematics and Statistics, and Department of Biomedical
Engineering, The University of Melbourne, Parkville, Victoria 3010, Australia
School of Medicine, Faculty of Medicine, Dentistry and Health Sciences, The
University of Melbourne, Parkville, Victoria 3010, Australia

###### Abstract

The maintenance of tissue and organ structures during dynamic homeostasis is
often not well understood. In order for a system to be stable, cell renewal,
cell migration and cell death must be finely balanced. Moreover, a tissue’s
shape must remain relatively unchanged. Simple epithelial tissues occur in
various structures throughout the body, such as the endothelium, mesothelium,
linings of the lungs, saliva and thyroid glands, and gastrointestinal tract.
Despite the prevalence of simple epithelial tissues, there are few models
which accurately describe how these tissues maintain a stable structure.

Here, we present a novel, 3D, deformable, multilayer, cell-centre model of a
simple epithelium. Cell movement is governed by the minimisation of a bending
potential across the epithelium, cell-cell adhesion, and viscous effects. We
show that the model is capable of maintaining a consistent tissue structure
while undergoing self renewal. We also demonstrate the model’s robustness
under tissue renewal, cell migration and cell removal. The model presented
here is a valuable advancement towards the modelling of tissues and organs
with complex and generalised structures.

##  1 Introduction

The body consists of four basic tissue types: epithelial, connective, muscular
and nervous [20]. The epithelium line the internal and external surfaces of
the body, as well as cavities, and some organs and glands. These epithelia are
the functioning components of the tissues, capable of performing various
tasks, including protection, secretion, absorption, excretion, filtration,
diffusion and sensation, depending on tissue location. Epithelial tissues have
a high potential for malignancies, in the form of cancers, as well as many
other diseases, such as asthma, cardiac disease, and many viral induced
diseases. Therefore, to better understand how disease develops, a clear
understanding towards how healthy epithelia are maintained in homeostasis is
first required.

Epithelia can be categorised in a number of ways, with one possible way being
the organisation, shape and function of the constituent cells. In terms of
cell shape, epithelial cells can be: squamous, cells with a width-to-height
ratio greater than 1; cuboidal, cells with a width-to-height ratio
approximately equal to 1; or columnar, cells with a width-to-height ratio less
than 1. Tissue organisation is categorised as either simple, if the epithelium
contains only a single layer of epithelial cells; stratified, if the
epithelium contains two or more layers of epithelial cells; or
pseudostratified, when the epithelium is made up of a single layer of columnar
cells of non uniform width along the long axis [20]. Simple epithelial tissues
are common throughout the human body, and occur in various geometric
structures, depending on the constituent cells’ shape. See Table 1 for a
schematic of simple epithelial tissue types, and some examples of their
locations within the body.

Epithelium type | Location | Reference
---|---|---
Simple squamous
|  Endothelium (capillary walls) Lung alveoli Mesothelium (peritoneum) |  [25] [7] [19]
Simple cuboidal
|  Ovary surface Renal tube lining Linings of the lungs Saliva glands Eyes Thyroid glands |  [21] [26] [18] [11] [15] [4]
Simple Columnar
|  Gastrointestinal tract Endocervix Fallopian tubes |  [24] [29] [27]
Table 1:  Schematic of simple epithelial tissues, and their locations within
the body.

Simple epithelia reside on a basement membrane, supported by the connected
stroma tissue [20]. They consist of self-renewing cells, which vary in cell
cycle duration depending on the tissue’s location within the body. Despite the
frequency of simple epithelial tissues throughout the body, the properties
leading to their particular geometric structures are not well understood. To
better understand their function, many different mathematical models of simple
epithelia exist. For example, models can incorporate cells as discrete
interacting agents. This is typically done on a 2D planar geometry [3], or a
fixed 3D geometry [12]. While valuable and insightful, these approaches
provide limited information toward the geometric structure of the epithelium.
Another approach is to model the epithelium as a continuum string (2D) [2] or
sheet (3D) [17] of cells. However, the insights gained into the geometric
structure from deformable continuum models comes with a loss of information
into cell topology.

In this work, we present a novel model of realistic 3D deformations of simple
epithelium. The aim of this model is to be able to accurately incorporate both
realistic tissue geometries and cell topology. The model is generalised, with
the goal of being customisable, based on the modelling applications.

##  2 Model of simple epithelium

Here, we present a mathematical model of a simple epithelium. We model the
tissue using a multicellular model, where cells are represented by their
centres, which are free to move in space [22]. We use multiple layers of
cells: an epithelial layer, in which cells are able to proliferate, and a
stromal layer, which supports the epithelium, in which cells are
differentiated and do not proliferate.

###  2.1 Biomechanical model

We use a lattice-free, cell-centred model, to describe the biomechanical
forces experienced by cells, with $\mathbf{r}_{i}$ denoting the position of
cell centre $i$ [23]. We denote the net force on a given cell $i$,
$\mathbf{F}_{i}^{\text{Net}}$. Here, we include forces due to neighbouring
interactions, $\mathbf{F}_{i}^{\text{Interaction}}$, viscous forces,
$\mathbf{F}_{i}^{\text{Viscous}}$ and lastly the bending force,
$\mathbf{F}_{i}^{\text{Bending}}$. The net force is therefore

$\displaystyle\mathbf{F}_{i}^{\text{Net}}=\mathbf{F}_{i}^{\text{Interaction}}+\mathbf{F}_{i}^{\text{Viscous}}+\mathbf{F}_{i}^{\text{Bending}},\quad\forall
i.$ (1)

We assume that cell motion is over-damped, due to the highly viscous cellular
environment [10]. Ignoring inertial terms, the net force on any given cell-
centre is zero, i.e. $\mathbf{F}_{i}^{\text{Net}}=\mathbf{0}$. We model
elastic neighbour interactions via a spring potential energy,
${U}^{\text{Interaction}}$. We first write the interaction potential for a
single cell $j$ as $P_{j}$:

$\displaystyle P_{j}=\sum_{n\in
N_{j}}\frac{1}{2}p_{jn}\left(s_{jn}-|\mathbf{r}_{jn}|\right)^{2},$ (2)

where $N_{j}$ is the set of the first-order neighbours, given by the Delaunay
triangulation, $0\leq p_{jn}\leq 1$ describes how readily cells $j$ and $n$
interact, $\mathbf{r}_{jn}$ is the displacement between cell-centres $j$ to
$n$, and $s_{jn}$ is the rest separation between cells $j$ and $n$ in the
absence of other external forces. The interaction potential of the tissue is

$\displaystyle{U}^{\text{Interaction}}=\sum_{\forall j}k_{j}P_{j},$ (3)

where $k_{j}\geq 0$ is the spring constant and denotes how readily cell $j$
adheres to its neighbours. The interaction force experienced by cell $i$ is
then defined as the minimisation of ${U}^{\text{Interaction}}$:

$\displaystyle\mathbf{F}_{i}^{\text{Interaction}}=-\nabla_{i}{U}^{\text{Interaction}},\quad\forall
i.$ (4)

We note that this is equivalent to a conventional spring force in the form
$\mathbf{F}_{i}=\sum_{n\in
N_{i}}\kappa_{in}\left(|\mathbf{r}_{in}|-s_{in}\right)\hat{\mathbf{r}}_{in}$,
where $\kappa_{in}$ is the spring constant between cells $i$ and $n$. The
viscous force simply opposes the direction of motion

$\displaystyle\mathbf{F}_{i}^{\text{Viscous}}=-\nu_{i}\frac{d\mathbf{r}_{i}}{dt},\quad\forall
i,$ (5)

where $\nu_{i}>0$ is the drag coefficient of cell $i$.

###  2.2 Bending force

If we want to examine the shape of a given tissue, we first require a measure
of “shape”. A common technique to quantify the shape of a triangulated surface
is the discrete Gaussian curvature [28], defined locally at each vertex $j$ in
the surface $S$ as:

$\displaystyle G_{j}$ $\displaystyle=\frac{2\pi-\phi_{j}}{A_{j}},\quad\forall
j,$ (6)

where $\phi_{j}=\sum_{n\in T_{j}}\theta_{n}$ is the angle sum at $j$, $T_{j}$
are the set of vertices that share an edge with vertex $j$, $\theta_{n}$ is
the angle of the $n^{\text{th}}$ triangular-element at vertex $j$, and
$A_{j}=\sum_{n\in T_{j}}a_{n}$ is the area contribution of the element that
contains vertex $j$ on the surface $S$, as shown by Figure 1(a) where $a_{n}$
is calculated as one third of the area of triangular-element $n$.

For a hexagonally-packed set of points, there are three possible
configurations for the element shown in Figure 1(a). We depict these
configurations in Figure 1(b)–1(d). If we consider the element shown in Figure
1(a) as a single element, we can say the element is elliptic if the discrete
Gaussian Curvature $G_{j}>0$, flat if $G_{j}=0$, and hyperbolic if $G_{j}<0$.

(a)

(b) $G_{j}>0$

(c) $G_{j}=0$

(d) $G_{j}<0$

Figure 1:  1(a) Local element used to calculate the discrete Gaussian
curvature. Possible element configurations, for different possible values of
$G_{j}$: an elliptical element shown in 1(b), a flat element shown in 1(c) and
hyperbolic element shown in 1(d).

Since our model defines cell connectivity within the tissue via a Delaunay
triangulation, we utilise this triangulation and consider only the epithelial
monolayer to defined a surface $S^{E}$. Using this surface, we classify the
global shape of the tissue through a bending potential, $U^{\text{Bending}}$,
as:

$\displaystyle U^{\text{Bending}}$ $\displaystyle=\sum_{\forall
j}\beta_{j}|G_{j}|^{\alpha},$ (7)

where $\alpha>0$ is the bending exponent and $\beta_{j}$ is the magnitude of
the bending potential at cell $j$ with:

$\displaystyle\beta_{j}=\begin{cases}\beta,&\text{if }j\in N^{E},\\\
0,&\text{otherwise},\end{cases}$ (8)

where $\beta\in\mathbb{R}$ and $N^{E}$ is the set of epithelial cells within
the tissue. We now minimise the bending potential, which will result in the
epithelial monolayer surface deforming to minimise the global curvature. In
this instance, minimising the bending potential results in a flat epithelial
monolayer surface, $S^{E}$. We therefore write the bending force on epithelial
cell centre $i$ as:

$\displaystyle\mathbf{F}^{\text{Bending}}_{i}$
$\displaystyle=-\nabla_{i}U^{\text{Bending}},\quad\forall i.$ (9)

For a detailed explanation of the exact form of the gradient of the bending
force, see SI.1. We can now write the equations of motion as:

$\displaystyle\nu_{i}\frac{d\mathbf{r}_{i}}{dt}=-\nabla_{i}\left(\sum_{\forall
j}k_{j}P_{j}+\sum_{\forall j}\beta_{j}|G_{j}|^{\alpha}\right),\quad\forall i.$
(10)

###  2.3 Tissue geometry

To best replicate in vitro and in vivo conditions, we simulate the tissue on a
periodic domain, to mimic the behaviour of a much larger tissue [16, 14, 1].
We impose periodic boundaries in the $x$ and $y$ axes. We allow the $z$ axis
to be a free boundary at both the top and bottom. Below, the tissue is
supported by the extracellular matrix, and above, the tissue is exposed to the
open lumen. However, since our simulation method uses a Delaunay triangulation
to define cell connectivity, we incorporate passive ghost nodes above and
bellow the tissue to prevent long edges forming between cells which cannot
possibly be connected to one another [23]. An example of the tissue
configuration with ghost nodes is shown in Figure 2.

The cell packing we use to initialise the tissue is called hexagonal close
packed. With a hexagonal close packed tissue, cells are initially separated by
their rest length, and therefore artefacts due to the tissue relaxing are not
present. We specify the tissue size as $\left[N_{x},N_{y},N_{z},N_{g}\right]$
where $N_{x}$ is the number of cells packed along the $x$ axis, $N_{y}$ the
number of cells packed along the $y$ axis, $N_{z}$ the number of cell layers,
and $N_{g}$ the number of ghost node layers to cap the tissue. $N_{x}$ and
$N_{y}$ give the spatial domain size as $\Omega_{x}=[0,N_{x}]$ and
$\Omega_{y}=[0,\frac{\sqrt{3}}{2}N_{y}]$ to ensure hexagonal packing of cells.
Since we always require that an epithelial monolayer is present, the number of
stromal cell layers is $N_{z}-1$. For example, the tissue in Figure 2 is of
size $\left[N_{x},N_{y},N_{z},N_{g}\right]=\left[8,10,2,1\right]$.

Figure 2:  Initial tissue set up of size
$\left[N_{x},N_{y},N_{z},N_{g}\right]=\left[8,10,2,1\right]$. Epithelial cells
are shown in pink, stromal cells in dark grey and ghost nodes are light
(transparent) grey.

###  2.4 Cell turnover

We choose proliferative cells to have a uniform cell cycle duration, $T\sim
U(11,13)$, with a mean of 12 hours [5, 8]. When a proliferative cell completes
its cell cycle, it is labelled a parent cell, with position $\mathbf{r}^{p}$.
The parent cell proliferates, producing two daughter cells, with positions
$\mathbf{r}_{i}$ and $\mathbf{r}_{j}$, separated by distance $\varepsilon>0$
apart. Following the division event, the parent cell is replaced by the two
daughter cells placed $\varepsilon$ apart on the plane defined by the
neighbours of the parent cell and the parent cell itself, spanned by the
orthonorma vectors $\mathbf{\hat{u}}$ and $\mathbf{\hat{v}}$. We can write the
daughter cell positions as:

$\displaystyle\mathbf{r}_{i}=\mathbf{r}^{p}+\frac{\varepsilon}{2}\left(\cos\theta\mathbf{\hat{u}}+\sin\theta\mathbf{\hat{v}}\right),\quad\mathbf{r}_{j}=\mathbf{r}^{p}-\frac{\varepsilon}{2}\left(\cos\theta\mathbf{\hat{u}}+\sin\theta\mathbf{\hat{v}}\right),$
(11)

where $\theta\in\left[0,2\pi\right)$ is drawn uniformly. We call this form of
division planar cell division. During a cell’s first hour of cell cycle, we
take the standard approach to increase the length separation between cell-
centres $i$ and $j$ uniformly as:

$\displaystyle
s_{ij}(t)=\begin{cases}\varepsilon+\tau_{i}\left(1-\varepsilon\right),\quad&\tau_{i}\leq
1,\text{ and $i,j$ share a parent cell},\\\
1,\quad&\text{otherwise},\end{cases}$ (12)

where $\tau_{i}\geq 0$ is the age of cell $i$ at time $t$ [22, 23].

We also know that if an epithelial cell detaches from the stroma, it is shed
tissue and dies, a process known as anoikis [30, 6, 13]. To account for this,
if an epithelial cell loses contact with the tissue, we remove it from the
simulation. Lastly, if cell $i$ is marked as being apoptotic, over the next
hour, the cell shrinks by changing the length separation as:

$\displaystyle
s_{ij}(t)=\begin{cases}1-2\tau^{A}_{i},\quad&0\leq\tau^{A}_{i}\leq\frac{1}{2}\\\
0,\quad&\frac{1}{2}<\tau^{A}_{i}\leq 1,\\\ 1,\quad&\text{if $j$ is a stromal
cell},\end{cases}$ (13)

where $\tau^{A}_{i}$ is the time since apoptosis began for cell $i$.

###  2.5 Implementation

To simulate the tissue, we choose $\alpha=1.01$ and $\beta=4$ for the bending
force parameters, as these provide the desired dynamics to maintain a flat,
epithelial monolayer surface, see SI.2. We also choose a spring constant,
$k_{i}=20$ and $p_{ij}=1$ if $i$ and $j$ are the same cell type, otherwise
$p_{ij}=0.5$, for the interaction model, and $\nu_{i}=2$ for cell viscosity.
We solve the equations of motion numerically using a Forward Euler method,
with a step size of $\Delta t=0.001$ hrs. Table 2 contains a summary of the
parameters used for the remainder of this paper. We implement the models used
within this paper in the Chaste environment, which is an open source C++
library, available at https://chaste.cs.ox.ac.uk [9]. The code developed for
this paper, is freely available at
https://github.com/DGermano8/ModelingRealistic3DDeformationsOfEpitheliumInDynamicHomeostasis.git.

Parameter | Value | Units | Reference
---|---|---|---
$\alpha$ | 1.01 | - | SI.2
$\beta$ | 4 | CD2 $\text{rad}^{-\alpha}$ | SI.2
$k_{i}$ | 20 | CM hrs-2 | [22]
$p_{ij}$ | $\begin{cases}1,\,\text{if $i,\,j$ are same cell type}\\\ 0.5,\,\text{otherwise}\end{cases}$ | - | -
$\nu_{i}$ | 2 | CM hrs-1 | [22]
$\Delta t$ | 0.001 | hrs | -
Table 2:  Parameters used for in silico tissue simulations. Cell diameter (CD)
is the natural spatial unit of the system and cell mass (CM) is the natural
mass units of the system.

An example of a tissue used initialise with size
$\left[N_{x},N_{y},N_{z},N_{g}\right]=\left[8,10,2,1\right]$ is shown in
Figure 2.

To measure how flat a given tissue is, we calculate the mean point curvature,
$\bar{G}$ as

$\displaystyle\bar{G}=\frac{1}{\left|N^{E}\right|}\sum_{j\in N^{E}}G_{j}.$
(14)

##  3 Results

We use our model to simulate three in silico experiments:

  * 1.

the relaxation of a deformable non-renewing differentiated tissue,

  * 2.

a deformable renewing tissue, to demonstrate the model’s robustness while
undergoing tissue renewal, in comparison to the base model,

  * 3.

cell migration in deformable renewing tissue, to demonstrate the model’s
robustness while undergoing tissue renewal, cell migration and cell removal,
in comparison to a traditional fixed geometry model.

###  3.1 Deformable non-renewing differentiated tissue

We first consider a tissue of differentiated cells relaxing. To do this, we
take a tissue with hexagonal packing, at equilibrium, and perturb the
$(x,y,z)$ positions by $(\zeta_{x},\zeta_{y},\zeta_{z})$, where $\zeta_{i}\sim
U(-\frac{1}{4},\frac{1}{4})$, for $i\in{x,y,z}$. An example of a perturbed
tissue is shown in Figure 3(a). The tissue is then free to relax, resulting in
a tissue similar to that shown in Figure 3(b). We take the average of 20
realisations of the in silico tissue simulations, and calculate a 95%
confidence interval. We see from the mean point curvature in Figure 3(c) that
the curvature starts at the maximum, as the initial epithelial monolayer
surface is not flat, and decreases as the bending force flattens the surface
$S^{E}$. This demonstrates the model’s ability to deform a non-flat tissue to
flat over time.

(a) $t=0$ hrs

(b) $t=2$ hrs

(c)

Figure 3:  In silico tissue simulations over 2 hours. Initial positions
perturbed from equilibrium. We only show the differentiated epithelial cells
(purple), the stromal cells (grey) and the epithelial monolayer surface
$S^{E}$ shown via triangulation. 3(a) is sample tissue, perturbed from
equilibrium, 3(b) is the final tissue. 3(c) shows the mean point curvature
over 20 sample tissue simulations (black line), with a 95% confidence interval
(red region). For an example video, see SI.3

###  3.2 Deformable renewing tissue

We next consider a renewing tissue, with cell proliferation and removal, with
results shown in Figure 4. To maintain a uniform tissue density and prevent
overcrowding, when a cell proliferates, we mark another epithelial cell as
apoptotic, chosen at random, unless an anoikis event has just occurred.

First, we simulate the tissue with no bending force, which we will call the
base model, and then again with identical initial conditions, this time with
the bending force, which we will call the deformable model, with examples of
the final state of the tissue model shown in Figures 4(a) and 4(b). We compute
$\bar{G}$ for both the base model and the deformable model. We show how the
mean point curvature varies as cells proliferate in Figure 4(c). We see that
the base model develops kinks and is no longer flat, whereas with the
deformable model, the tissue remains flat. We can also see that the mean point
curvature for the two different tissues varies by up to 3 orders of magnitude,
with the deformable model having a significantly lower mean point curvature.
Lastly, we note that the deformable model contains distinct spikes in the mean
point curvature, which are not sustained, but disappear as the tissue
flattens. These spikes are caused by cells proliferating. We note that the
same spikes are also present in the base model, but appear less significant
due to the log scale. These results show that the deformable models ability to
maintain structure is robust while the tissue undergoes self renewal.

(a) Base model at $t=48$ hrs

(b) Deformable model at $t=48$ hrs

(c)

Figure 4:  Tissue simulations with proliferation and removal over 48 hours.
For simplicity, we only show the epithelial cells, the epithelial monolayer
surface $S^{E}$, and stromal cells. 4(a) and 4(b) are the final states of the
tissues for the base model and the deformable model respectively. 4(c) shows
how the mean point curvature of the tissues varies with time, as cells
proliferate. Blue shows the tissue for the base model, and purple the tissue
with the deformable model, along with 95% confidence intervals (red regions,
calculated over 5 unique simulations). For an example video, see SI.4.

###  3.3 Cell migration in deformable renewing tissue

We now consider both proliferative epithelial cells and differentiated
epithelial cells, to observe how our model behaves with cell migration. We
implement a simple position based model, where a cell’s position determines
its cell type. In the proliferative region, with
$\mathbf{r}_{i}\in\Omega_{P}$, the cell experiences high external signalling
factor and the cell is proliferative, otherwise a cell is differentiated. We
also implement a density dependent death model where, if a cell’s area
contribution is below a certain threshold (i.e. $A_{j}<A_{\text{crit}}$), and
the cell is located within the death region
$\mathbf{r}_{i}\in\Omega_{D}=\Omega^{L}_{D}\cup\Omega^{R}_{D}$, the cell is
marked as apoptotic and undergoes apoptosis over the next hour. Figure 5 shows
a diagram of the proliferative, migratory and death regions.

Figure 5:  Domain structure with the proliferative, $\Omega_{P}$, migratory, $\Omega_{M}=\Omega^{L}_{M}\cup\Omega^{R}_{M}$, and death, $\Omega_{D}=\Omega^{L}_{D}\cup\Omega^{R}_{D}$, regions. Parameter | Value | Units
---|---|---
$A_{\text{crit}}$ | $0.98\,A_{\text{Relaxed}}$ | CD2
$A_{\text{Relaxed}}$ | $\frac{\sqrt{3}}{2}$ | CD2
$\Omega_{P}$ | $\\{(x,y,z)\in\mathbb{R}^{3}:x\in[-1.5,1.5]\\}$ | (CD,CD,CD)
$\Omega_{M}$ | $\\{(x,y,z)\in\mathbb{R}^{3}:x\in[-5.5,-1.5)\cup(1.5,5.5]\\}$ | (CD,CD,CD)
$\Omega_{D}$ | $\\{(x,y,z)\in\mathbb{R}^{3}:x\in[-6,-5.5)\cup(5.5,6]\\}$ | (CD,CD,CD)
Table 3:  Parameters used for in silico tissue simulations. Cell diameter (CD)
is the natural spatial unit of the system. $A_{\text{Relaxed}}$ is the area of
a hexagon, which is the resting equilibrium of cells in hexagonal close
packing.

Using the domain structure shown in Figure 5, and the domain sizes and
critical area for density-dependent death given in Table 3, we simulate a
tissue of size $\left[N_{x},N_{y},N_{z},N_{g}\right]=\left[12,14,3,1\right]$
for 480 hours.

(a)

(b)

(c)

Count

Figure 6: An example of the tissue simulation at steady state. Figure 6(a)
shows hyper stimulated proliferative transient cells, shown in pink, migratory
differentiated cells, shown in purple, and dying cells, shown in red. Stromal
cells are shown in grey, and ghost nodes are not shown. Figures 6(b) and 6(c)
show epithelial cell distributions throughout the tissue at homeostasis. The
sample distributions are collated by sampling the tissue every hour for the
last 360 hours of the simulation to ensure homeostasis. For an example video,
see SI.5

To ensure the tissue is in dynamic homeostasis, we discard the first 120
hours, and consider the tissue from time $t=120hrs$ to $t=480hrs$. In Figures
6(b) and 6(c) we show the cell distribution for the tissue in homeostasis,
from 360 hourly samples of a simulation. In Figure 6(b) we show that the
tissue undergoes deformation in the vertical ($z$) axis with respect to the
migration and removal ($x$) axis. However, this deformation is small at
approximately 1CD across the 12CD migration and removal ($x$) axis. We can
also see that, with respect to the migration and removal ($x$) axis, cells are
uniformly distributed throughout. If we consider Figure 6(c) which shows the
cell distribution with the horizontal ($y$) and the vertical ($z$) axis. Here
we can see that cells occupy distinct strips within the tissue, which
correspond to the migration paths of the cells through the tissue. We note
that these distinct bands are a result of the domain specification chosen.

To observe how the model presented here behaves in comparison to previous
models of cell migration, we consider a tissue evolving on a fixed geometry
and compared it to the above simulation. The fixed geometry is such that cells
are free to move along the migration and removal axis, as well as the
horizontal axis, but epithelial cells are limited to a fixed vertical height,
$z_{F}$. Stromal cells and ghost nodes do not have any restrictions. We note
that this is similar to considering a purely 2D model of cell migration.
However, a 2D model is not sufficient as it does not account for the cell-cell
interactions between epithelial and stromal cells, which has the effect of an
increased effective viscosity in 2D. As such, we compare our model to a fixed
geometry model rather than a purely 2D model.

(a)

(b)

Figure 7: Distribution of migration times for fixed (grey) and deformable
(blue) models are shown in 7(a). The migration speeds along the cell migration
and removal ($x$) axis are shown in 7(b), for both the fixed (grey) and
deformable (blue) models, each with 95% confidence intervals. The background
colours of Figure 7(b) show the region the sample point is in.

Figure 7(a) shows the distribution of migration times for both the fixed and
deformable models. The migration time is defined as the time taken for a cell
to migrate throughout the migratory region, $\Omega_{M}$. We see that while
the deformable model has a longer tail, and is slightly skewed to lower
migration time, both tissues behave in a similar manner.

If we discretize the migration and removal axis, and look at a cells migration
speed within each region, we observe how cells are migrating throughout the
tissue, as shown in Figure 7(b), for both fixed (grey) and deformable (blue)
models. We can see that, within the proliferative region, $\Omega_{P}$, cell
migration speeds overlap. At the interface between the proliferative and
migratory region, the fixed model has a higher cell migration speed. However
throughout the remainder of the tissue, migrations speeds agree, within the
margin of error of the 95% confidence intervals.

The analysis shown in Figure 7 indicates that the deformable model does not
unduly influence the dynamics found in fixed geometry models of dynamic
homeostasis, while still allowing deformations to be included in such systems.
These results show that the deformable models ability to maintain structure is
robust while the tissue undergoes renewal, cell migration and cell removal.

##  4 Discussion

We have presented a novel model of a 3D deformable, dynamic tissue. Using a
cell-centre based Delaunay triangulation model, we described cell interactions
via a spring potential, and using the same triangulation, we developed a
bending potential to describe the shape of a tissue. Through the minimisation
of this bending potential, we displayed how we are able to maintain a flat
epithelial monolayer, undergoing proliferation and death, which resides on a
stromal tissue.

We then compared how a simulated tissue behaves both with and without a
bending force acting on the tissue. We found that when we do not include a
bending force, the tissue develops kinks and no longer remains flat. In
comparison, when we included the bending force, we were able to maintain a
flat epithelial monolayer, with a bending potential of up to 3 orders of
magnitude smaller than the tissue without a bending force.

Finally, we implemented the model on a tissue which experiences localised cell
proliferation, and a local density-dependent death model. We found that with
our model presented here, we are able to achieve cell migration while still
maintaining a relatively flat epithelial monolayer. We also compared our model
to a conventional 3D fixed geometry model, and found that the migration time
and migration speed follow similar distributions. From this, we can conclude
that the model presented here is well suited to capture the compartmental and
deformable nature of a realistic tissue in dynamic homeostasis.

Future avenues of study include incorporating a realistic curved tissues, and
to determine the conditions under which does the system remains stable. We
will then apply the model to some commonly studied tissues and organs. We
would also like to compare the results here with those obtained from continuum
models of tissue deformation and cell-based models of deformable tissues, to
draw similarities and differences.

## Acknowledgements

This research was supported by an Australian Government Research Training
Program (RTP) Scholarship (awarded to DPJG). S.T.J. is supported by the
Australian Research Council (project no. DE200100998).

## References

  * MIL [2021] Maintaining the proliferative cell niche in multicellular models of epithelia.  _Journal of Theoretical Biology_ , 527:110807, 2021.  ISSN 0022-5193.  doi: https://doi.org/10.1016/j.jtbi.2021.110807.
  * Almet et al. [2021] A. A. Almet, H. M. Byrne, P. K. Maini, and D. E. Moulton.  The role of mechanics in the growth and homeostasis of the intestinal crypt.  _Biomechanics and Modeling in Mechanobiology_ , 20(2):585–608, 2021.
  * An [2008] G. An.  Introduction of an agent-based multi-scale modular architecture for dynamic knowledge representation of acute inflammation.  _Theoretical Biology and Medical Modelling_ , 5(1):1–20, 2008.
  * Balasubramanian [2020] S. P. Balasubramanian.  Anatomy of the thyroid, parathyroid, pituitary and adrenal glands.  _Surgery (Oxford)_ , 38(12):758–762, 2020.
  * Barker et al. [2009] N. Barker, R. A. Ridgway, J. H. Van Es, M. Van De Wetering, H. Begthel, M. Van Den Born, E. Danenberg, A. R. Clarke, O. J. Sansom, and H. Clevers.  Crypt stem cells as the cells-of-origin of intestinal cancer.  _Nature_ , 457(7229):608–611, 2009.
  * Battini et al. [2006] L. Battini, E. Fedorova, S. Macip, X. Li, P. D. Wilson, and G. L. Gusella.  Stable knockdown of polycystin-1 confers integrin-$\alpha$2$\beta$1–mediated anoikis resistance.  _Journal of the American Society of Nephrology_ , 17(11):3049–3058, 2006.
  * Bonastre et al. [2016] E. Bonastre, E. Brambilla, and M. Sanchez-Cespedes.  Cell adhesion and polarity in squamous cell carcinoma of the lung.  _The Journal of Pathology_ , 238(5):606–616, 2016\.
  * Cheng and Leblond [1974] H. Cheng and C. Leblond.  Origin, differentiation and renewal of the four main epithelial cell types in the mouse small intestine v. unitarian theory of the origin of the four epithelial cell types.  _American Journal of Anatomy_ , 141(4):537–561, 1974.
  * Cooper et al. [2020] F. R. Cooper, R. E. Baker, M. O. Bernabeu, R. Bordas, L. Bowler, A. Bueno-Orovio, H. M. Byrne, V. Carapella, L. Cardone-Noott, J. Cooper, et al.  Chaste: cancer, heart and soft tissue environment.  _Journal of Open Source Software_ , 2020.
  * Dallon and Othmer [2004] J. C. Dallon and H. G. Othmer.  How cellular movement determines the collective force generated by the dictyostelium discoideum slug.  _Journal of Theoretical Biology_ , 231(2):203–222, 2004.
  * de Paula et al. [2017] F. de Paula, T. H. N. Teshima, R. Hsieh, M. M. Souza, M. M. S. Nico, and S. V. Lourenco.  Overview of human salivary glands: highlights of morphology and developing processes.  _The Anatomical Record_ , 300(7):1180–1188, 2017\.
  * Dunn et al. [2013] S.-J. Dunn, I. S. Näthke, and J. M. Osborne.  Computational models reveal a passive mechanism for cell migration in the crypt.  _PLoS One_ , 8(11), 2013.
  * Eisenhoffer et al. [2012] G. T. Eisenhoffer, P. D. Loftus, M. Yoshigi, H. Otsuna, C.-B. Chien, P. A. Morcos, and J. Rosenblatt.  Crowding induces live cell extrusion to maintain homeostatic cell numbers in epithelia.  _Nature_ , 484(7395):546–549, 2012.
  * Fletcher et al. [2013] A. G. Fletcher, J. M. Osborne, P. K. Maini, and D. J. Gavaghan.  Implementing vertex dynamics models of cell populations in biology within a consistent computational framework.  _Progress in Biophysics and Molecular Biology_ , 113(2):299–326, 2013.
  * Frost et al. [2014] L. S. Frost, C. H. Mitchell, and K. Boesze-Battaglia.  Autophagy in the eye: implications for ocular cell health.  _Experimental Eye Research_ , 124:56–66, 2014.
  * Germano and Osborne [2020] D. P. Germano and J. M. Osborne.  A mathematical model of cell fate selection on a dynamic tissue.  _Journal of Theoretical Biology_ , 2020.  doi: doi.org/10.1016/j.jtbi.2020.110535.
  * Hannezo et al. [2011] E. Hannezo, J. Prost, and J.-F. Joanny.  Instabilities of monolayered epithelia: shape and structure of villi and crypts.  _Physical Review Letters_ , 107(7):078104, 2011\.
  * Hermans and Bernard [1999] C. Hermans and A. Bernard.  Lung epithelium–specific proteins: characteristics and potential applications as markers.  _American Journal of Respiratory and Critical Care Medicine_ , 159(2):646–678, 1999.
  * Hiriart et al. [2019] E. Hiriart, R. Deepe, and A. Wessels.  Mesothelium and malignant mesothelioma.  _Journal of Developmental Biology_ , 7(2):7, 2019\.
  * Kahle and Frotscher [1976] W. Kahle and M. Frotscher.  _Color Atlas Textbook of Human Anatomy. Vol. 3. Nervous System and Sensory Organs_.  Thieme, 1976.
  * Katabuchi and Okamura [2003] H. Katabuchi and H. Okamura.  Cell biology of human ovarian surface epithelial cells and ovarian carcinogenesis.  _Medical Electron Microscopy_ , 36(2):74–86, 2003\.
  * Meineke et al. [2001] F. A. Meineke, C. S. Potten, and M. Loeffler.  Cell migration and organization in the intestinal crypt using a lattice-free model.  _Cell Proliferation_ , 34(4):253–266, 2001.
  * Osborne et al. [2017] J. M. Osborne, A. G. Fletcher, J. M. Pitt-Francis, P. K. Maini, and D. J. Gavaghan.  Comparing individual-based approaches to modelling the self-organization of multicellular tissues.  _PLoS computational biology_ , 13(2):e1005387, 2017.
  * Reed and Wickham [2009] K. K. Reed and R. Wickham.  Review of the gastrointestinal tract: from macro to micro.  In _Seminars in Oncology Nursing_ , volume 25, pages 3–14. Elsevier, 2009.
  * Stolz and Sims-Lucas [2015] D. B. Stolz and S. Sims-Lucas.  Unwrapping the origins and roles of the renal endothelium.  _Pediatric Nephrology_ , 30(6):865–872, 2015\.
  * Thakur and Tiwari [2019] R. Thakur and A. Tiwari.  Determination of the time since death by histological changes in distal convoluted tubule in human kidneys.  _Int J Anat Res_ , 7(4.2):7086–91, 2019.
  * Wira et al. [2005] C. R. Wira, K. S. Grant-Tschudy, and M. A. Crane-Godreau.  Epithelial cells in the female reproductive tract: a central role as sentinels of immune protection.  _American Journal of Reproductive Immunology_ , 53(2):65–76, 2005.
  * Xu and Xu [2009] Z. Xu and G. Xu.  Discrete schemes for gaussian curvature and their convergence.  _Computers & Mathematics with Applications_, 57(7):1187–1195, 2009.
  * Yi et al. [2013] T. J. Yi, B. Shannon, J. Prodger, L. McKinnon, and R. Kaul.  Genital immunology and hiv susceptibility in young women.  _American Journal of Reproductive Immunology_ , 69:74–79, 2013.
  * Yin et al. [2022] J. Yin, J. Wang, X. Zhang, Y. Liao, W. Luo, S. Wang, J. Ding, J. Huang, M. Chen, W. Wang, et al.  A missing piece of the puzzle in pulmonary fibrosis: anoikis resistance promotes fibroblast activation.  _Cell & Bioscience_, 12(1):1–19, 2022.

## Supplementary Information

##  SI.1 Equations of motion

To calculate the exact form of the bending potential, we take the gradient
with respect to cell $i$ and note that since only the triangular components of
$\phi_{j}$ that contain cell centre $i$ will contribute:

$\displaystyle\mathbf{F}^{\text{Bending}}_{i}$
$\displaystyle=-\nabla_{i}\sum_{j\in
N^{E}}\beta_{j}\left|\frac{2\pi-\phi_{j}}{A_{j}}\right|^{\alpha},$
$\displaystyle\forall i\in N^{E},$ (SI.1) $\displaystyle=-\sum_{j\in
N^{\text{E}}_{i}}\beta_{i}\nabla_{i}\left|\frac{2\pi-\phi_{j}}{A_{j}}\right|^{\alpha},$
$\displaystyle\forall i\in N^{E},$ (SI.2) $\displaystyle=-\alpha\sum_{j\in
N^{\text{E}}_{i}}\beta_{i}\text{sgn}\left(\frac{2\pi-\phi_{j}}{A_{j}}\right)\left|\frac{2\pi-\phi_{j}}{A_{j}}\right|^{\alpha-1}\,\left[\frac{A_{j}\nabla_{i}\phi_{j}-(2\pi-\phi_{j})\nabla_{i}A_{j}}{A_{j}^{2}}\right],$
$\displaystyle\forall i\in N^{E},$ (SI.3)

where $N^{\text{E}}_{i}$ are the first-order, epithelial neighbours of cell-
centre $i$. From here, we see how $\beta_{i}$ describes the magnitude of the
bending force, and $\alpha$ the sensitivity of the bending force, which are
further discussed in SI.2. To maintain our epithelial monolayer surface,
$S^{E}$, the bending force is only applied to the epithelial cells. The
gradient of $\phi_{j}$, with respects to cell $i$ is:

$\displaystyle\mathbf{\nabla}_{i}\phi_{j}=\begin{cases}\frac{-1}{\sqrt{1-\cos^{2}(\theta_{mjn})}}\left[\left(\frac{\cos(\theta_{mjn})}{|\mathbf{r}_{jm}|}-\frac{1}{|\mathbf{r}_{jn}|}\right)\mathbf{\hat{r}}_{jm}+\left(\frac{\cos(\theta_{mjn})}{|\mathbf{r}_{jn}|}-\frac{1}{|\mathbf{r}_{jm}|}\right)\mathbf{\hat{r}}_{jn}\right],\quad&\text{if
}j=i,\\\
\frac{-1}{\sqrt{1-\cos^{2}(\theta_{jim})}}\frac{1}{|\mathbf{r}_{ij}|}\left[\mathbf{\hat{r}}_{im}-\mathbf{\hat{r}}_{ij}\cos(\theta_{jim})\right],\quad&\text{if
}j\in N^{E}_{i},\\\ \mathbf{0},\quad&\text{otherwise},\end{cases}$ (SI.4)

where $\theta_{mjn}$ is the angle between $\mathbf{r}_{jm}$ and
$\mathbf{r}_{jn}$. Similarly, the gradient of $A_{j}$, with respect to cell
$i$ is:

$\displaystyle\mathbf{\nabla}_{i}A_{j}=\begin{cases}\frac{1}{4A_{j}}\left[\psi_{n}\left(\mathbf{r}_{nm},\mathbf{r}_{mj}\right)+\psi_{m}\left(\mathbf{r}_{mn},\mathbf{r}_{nj}\right)+\psi_{j}\left(\mathbf{r}_{mn},\mathbf{r}_{mn}\right)\right],\quad&\text{if
}j\in N^{E}_{i},\\\ \mathbf{0},\quad&\text{otherwise},\end{cases}$ (SI.5)

where
$\psi_{l}\left(\mathbf{r}_{ab},\mathbf{r}_{cd}\right)=\mathbf{r}_{l}\left(\mathbf{r}_{ab}\cdot\mathbf{r}_{cd}\right)-\mathbf{r}_{l}\cdot\left(\mathbf{r}_{ab}\circ\mathbf{r}_{cd}\right)$,
and $\mathbf{r}_{ab}\circ\mathbf{r}_{cd}$ is the Hadamard product between
vectors $\mathbf{r}_{ab}$ and $\mathbf{r}_{cd}$. From here, we can rearrange
to give the following equations of motion:

$\displaystyle\nu_{i}\frac{d\mathbf{r}_{i}}{dt}=\sum_{j\in
N_{i}}k_{i}p_{ij}\left(|\mathbf{r}_{ij}|-s_{ij}\right)\mathbf{\hat{r}}_{ij}-\sum_{j\in
N^{\text{E}}_{i}}\alpha\beta_{i}\text{sgn}\left(\frac{2\pi-\phi_{j}}{A_{j}}\right)\left|\frac{2\pi-\phi_{j}}{A_{j}}\right|^{\alpha-1}\,\left[\frac{A_{j}\nabla_{i}\phi_{j}-(2\pi-\phi_{j})\nabla_{i}A_{j}}{A_{j}^{2}}\right],\quad\forall
i.$ (SI.6)

We now define the magnitude of the bending potential at cell $j$, $\beta_{j}$,
in terms of the spring constant $k_{j}$, and a proportionality constant,
$\rho_{j}$, which shares the same non-zero properties as $\beta_{j}$, and
write:

$\displaystyle\Rightarrow\nu_{i}\frac{d\mathbf{r}_{i}}{dt}$
$\displaystyle=\sum_{j\in
N_{i}}k_{i}\left\\{p_{ij}\left(|\mathbf{r}_{ij}|-s_{ij}\right)\mathbf{\hat{r}}_{ij}-\alpha\rho_{i}\,\text{sgn}\left(\frac{2\pi-\phi_{j}}{A_{j}}\right)\left|\frac{2\pi-\phi_{j}}{A_{j}}\right|^{\alpha-1}\,\left[\frac{A_{j}\nabla_{i}\phi_{j}-(2\pi-\phi_{j})\nabla_{i}A_{j}}{A_{j}^{2}}\right]\right\\},\quad\forall
i.$ (SI.7)

##  SI.2 Bending force calibration

To calibrate the bending force in our model, we consider a test element at
equilibrium, and perturb the $z$-position of the centre cell by $\delta z$. We
also consider the effect compression has on the model, where $C_{x}$ and
$C_{y}$ are the axis scalings in the $x$ and $y$ components respectively, and
where $C_{x}=C_{y}=1$ is no compression. The test element is shown in Figure
SI.1. We see that since planar forces in $\mathbf{F}_{i}^{\text{Interaction}}$
and $\mathbf{F}_{i}^{\text{Bending}}$ each individually balance, the
contribution due to these forces will be in the $z$-component only.
Specifically $\mathbf{F}_{i}^{\text{Bending}}$ will act downwards, and
$\mathbf{F}_{i}^{\text{Interaction}}$ will act upwards. We therefore only need
to consider the magnitudes of $\mathbf{F}_{i}^{\text{Bending}}$ and
$\mathbf{F}_{i}^{\text{Interaction}}$ .

Figure SI.1:  Test element used to calibrate model. $C_{x}$ and $C_{y}$ are
the amount of compression in the $x$ and $y$ components respectively, and
$\delta z$ the $z$ displacement from equilibrium.

Before we perform our calibration, we scale the bending force with respect to
the interaction force, writing $\beta_{i}$ in terms of $k_{i}$ and a
proportionality variable, $\rho_{i}$, which shares the same non-zero
properties as $\beta_{i}$, and define $\beta=k_{i}\rho$. We can then rewrite
Equation (10) as:

$\displaystyle\nu_{i}\frac{d\mathbf{r}_{i}}{dt}=-\nabla_{i}\sum_{\forall
j}k_{j}\left(P_{j}+\rho_{j}|G_{j}|^{\alpha}\right),\quad\forall i.$ (SI.8)

We now use the test element in Figure SI.1 to determine the exponent
parameter, $\alpha$, and the proportionality constant, $\rho$. Since we will
have equal compression in the $x$ and $y$ axes, we denote $C=C_{x}=C_{y}$.

We specify $k_{i}=20$ and $p_{ij}=1$, which gives the interaction force shown
in Figure 2(a), and the bending force for $\rho=0.2$ and $\alpha=1.01$ shown
in Figure 2(b). We can find the sum
$|\mathbf{F}^{\text{Interaction}}|-|\mathbf{F}^{\text{Bending}}|$, as shown in
Figure SI.3, and choose $\alpha$ and $\rho$ such that we achieve the desired
behaviour of maintaining a flat epithelial monolayer of cells.

$|\mathbf{F}^{\text{Interaction}}|$

(a)

$-|\mathbf{F}^{\text{Bending}}|$

(b)

Figure SI.2:  Figure 2(a) shows the cell-cell interaction force on the centre
cell of the test element, with $k_{i}=20$ and $p_{ij}=1$. Figure 2(b) shows
the bending force on the same cell with $\rho=0.2$ and $\alpha=1.01$.

Variation in $|\mathbf{F}^{\text{Interaction}}|-|\mathbf{F}^{\text{Bending}}|$
with $\rho$ and $\alpha$

$\rho=0.2$

$\alpha=1.01$

(a)

$\alpha=1.25$

(b)

$\alpha=1.5$

(c)

$\rho=0.1$

(d)

(e)

(f)

Figure SI.3:
$|\mathbf{F}^{\text{Interaction}}|-|\mathbf{F}^{\text{Bending}}|$ for varying
values of $\rho$ and $\alpha$. For $\alpha=1.01$, 3(d) and 3(a) show the
element comes to a near flat configuration. For $\alpha=1.25$, 3(e) and 3(b)
show the element will flatten out, but still come to rest at a height above
flat. For $\alpha=1.5$, 3(f) and 3(c) show that as the element is compressed,
the cell will always sit above flat, at a height dictated by the zero contour.

Figure SI.3 shows that as we decrease $\alpha$, we control the nullcline which
the element will achieve, if we were to move the centre cell only. We need the
smallest possible $\alpha$ to maintain a flat layer when the tissue is
compressed. However, if we choose $\alpha$ too small, the stability of the
numerical solver is compromised. This can be managed by choosing the
appropriate $\rho$ to ensure numerical stability is maintained. For this
paper, we specify $\rho=0.2$ and $\alpha=1.01$.

##  SI.3 Video: deformable non-renewing differentiated tissue example

The is a simulation sample showing the relaxation of a deformable non-renewing
differentiated tissue. Visualised are the epithelial cells (purple -
differentiated), the stromal cells, and the triangulation the bending force
acts upon.
https://drive.google.com/file/d/1dOSyDGMHbestnN6jFnaYomwx4XfH0497/view?usp=sharing.

##  SI.4 Video: deformable renewing tissue example

The is a simulation sample of a deformable renewing tissue. It shows that as
the tissue renews, the epithelial monolayer is maintained. Visualised are the
epithelial cells (pink - proliferative and red - apoptotic), the stromal
cells, and the triangulation the bending force acts upon.
https://drive.google.com/file/d/1aV53ALxFsaygoVmHa4t6WfSULZCl8rcs/view?usp=sharing.

##  SI.5 Video: cell migration in deformable renewing tissue example

The is a simulation sample of cell migration in a deformable renewing tissue.
It shows that as the cells proliferate in the proliferative region, they then
migrate and become differentiated, eventually undergoing apoptosis when they
become too compressed. The simulation also demostrates how, even under this
process, the tissues structure is maintained. Visualised are the epithelial
cells (pink - proliferative, purple - differentiated and red - apoptotic), the
stromal cells, and the triangulation the bending force acts upon.
https://drive.google.com/file/d/1ElUK5eCSEkcpuDflnI3UVBMKaLFLyA5V/view?usp=sharing.