# Solving ODEs vs VoronoiMesh

I'm using VoronoiMesh to model cells in an epithelium. The concentration of certain proteins within the cells are determined by solving a system of ODEs per cell. My goal is to model emerging patterns based on the concentrations of such proteins in each cell, where the colour of each cell is determined by the concentration of a specific protein.

I have two questions regarding the implementation of this simulation:

1. Is it possible to change the colour of each cell in a Voronoi mesh, given, for example, an array of the solutions of such concentrations? I'm aware of the MeshPrimitives command, but I'm not sure how to efficiently use it.

2. Also, these concentrations depend on the protein concentrations of neighbouring cells (these are cell-cell signalling dynamics). Is it possible to detect neighbouring cells so that I can incorporate their protein concentrations in the ODE system for a specific cell? Say, something like NeighbourQ? Maybe changing the whole thing into a graph could help, but I'd appreciate any suggestions.

Any ideas? Thank you.

The following example makes use of my package, IGraph/M.

Needs["IGraphM"]


First we generate some random points to use as Voronoi centres. We apply a few steps of Lloyd relaxation to make them more equi-distant.

pts = RandomReal[{-1, 1}, {100, 2}];
Do[
pts = PropertyValue[{VoronoiMesh[pts, {{-1, 1}, {-1, 1}}], {2, All}}, MeshCellCentroid],
{2}
]


This is our mesh, representing the cells:

mesh = VoronoiMesh[pts, {{-1, 1}, {-1, 1}},
MeshCellStyle -> {1 -> White}]


We are going to simulate diffusion between the cells. For this, we need the Laplacian matrix (also called Kirchhoff matrix) of their neighbourhood graph.

ag = IGMeshCellAdjacencyGraph[mesh, 2];
km = IGKirchhoffMatrix[ag]


IGraph/M provides IGMeshCellAdjacencyGraph to build the neighbourhood graph of mesh cells of various dimensions, IGMeshCellAdjacencyMatrix to build the corresponding adjacency matrix (faster), and IGKirchhoffMatrix to get the graph Laplacian. There is more information about these functions in the package's documentation (IGDocumentation[]).

Let us create random initial concentrations:

initC = RandomReal[{0, 1}, MeshCellCount[mesh, 2]];


Then solve the diffusion equation using NDSolve:

diffConst = 1.0;
solfun = NDSolveValue[{c'[t] == - diffConst * km.c[t], c[0] == initC}, c, {t, 0, 1}]


We can animate the solution like this:

Animate[
SetProperty[{mesh, {2, All}}, MeshCellStyle -> ColorData["Rainbow"] /@ solfun[t]],{t, 0, 1}
]


• Thank you. This looks promising, I will take a deeper look. Just one quick question: my ultimate goal is to deploy an interactive interface (using Manipulate) to the Cloud and be able to share it with other researchers, most of which will only install the WolframPlayer to use it. How compatible is this package with this player? Do people need to install it manually on Mathematica (this is something I'd like to avoid)? I wonder if commands like Get or Need work inside a Manipulate. Dec 28 '19 at 15:45
• @samwolfe It does not work with the free version of the Wolfram Player. It requires a full version of Mathematica. I have no experience with the paid version of the Wolfram Player. It does not work in the cloud because they do not allow arbitrary binaries there. Installing it into a full version of Mathematica is quite trivial though (see the instructions). Dec 28 '19 at 16:05
• I will let someone else post an answer that does not depend on my package. Dec 28 '19 at 16:07
• Thanks! One more thing: instead of MeshCellStyle -> ColorData["Rainbow"]/@ solfun[t], I wish the solution determines the opacity of a certain RGB colour, that is, something like RGBColor[0.5, 0.34, 0.5,solfun[t]]`, which doesn't work. How can I do that? Dec 28 '19 at 19:22
• Also, how do you guarantee that the right cell is coloured? It seems I'm getting the colouring in the wrong cells. Dec 28 '19 at 20:50