Imagine I want to model 2D cells moving in a tissue. Consider, as an example, Lloyd's relaxation algorithm. Using VoronoiMesh I'm able to mimic some cell movement. Now, if I colour one of the cells in a different manner, I can track it as the mesh moves.
However, this doesn't work as expected, since VoronoiMesh yields different cell ordering for different sets of points and therefore I get the wrong colouring in the wrong cell.
Which was obtained by the following code, where the function VorR is inspired by this answer
rel = Function[{pts, zmp},
Block[{cells},
cells = MeshPrimitives[
VoronoiMesh[pts, {{-zmp, zmp}, {-zmp, zmp}}], "Faces"];
RegionCentroid /@
cells[[SparseArray[Outer[#2@#1 &, pts, RegionMember /@ cells, 1],
Automatic, False]["NonzeroPositions"][[All, 2]]]]]];
VorR = Function[{pt, s}, Module[{pts2, vor, vcells, mesh},
pts2 =
Flatten[Table[
TranslationTransform[{ 2 s i, 2 s j}][pt], {i, -1, 1}, {j, -1,
1}], 2];
vor = VoronoiMesh[pts2];
vcells = Catenate[NearestMeshCells[{vor, 2}, #] & /@ pt];
mesh = MeshRegion[MeshCoordinates[vor], MeshCells[vor, vcells]]
]];
n = 20; T = 50;
pts = RandomReal[{-1, 1}, {n, 2}];
val = ReplacePart[ConstantArray[1, n], RandomInteger[n] -> 0.2];
vl = {pts};
vorl = {VorR[vl[[1]], 1]};
For[i = 1, i <= T, i++,
sca = (i 0.01 + 1);
vl = Append[vl, rel[Last[vl], 1.05]];
vorl = Append[vorl, VorR[Last[vl]*sca, sca]];
];
colr = ConstantArray[Table[val[[j]], {j, n}], T];
cels = Table[
Table[Graphics[{RGBColor[0.5, 0.65, 0.5, colr[[j, p]]],
MeshPrimitives[vorl[[j]], 2][[p]]}], {p, n}], {j, T}];
ListAnimate[Table[Show[cels[[j]]], {j, T}]]
As it can be seen, the selected cell jumps between different places and this is due to the reordering of the cells (polygons) in MeshPrimitives[vorl[[i]],2], which occurs every time VoronoiMesh is applied. How can I solve this and get the correct ordering each time?
My attempt: Initially, I thought I could track the seeds updating positions (by Lloyd's algorithm) and associate each with a cell polygon, but that doesn't seem to work because the cells are not ordered by their generating seeds. Since movement is "relatively" slow, I then thought about tracking the position of the cell centroids in order to test, at each step, which cell permutation occurred and therefore correct the order at every step. And this seems to work, as seen by the following code (apologises for the "messy" code, there are definitely neater ways of writing this using Select for instance)
n = 20; T = 50;
pts = RandomReal[{-1, 1}, {n, 2}];
val = {ReplacePart[ConstantArray[1, n], RandomInteger[n] -> 0.2]};
vl = {pts};
vorl = {VorR[vl[[1]], 1]};
vlc = {Map[RegionCentroid, MeshPrimitives[vorl[[1]], 2]]};
val1 = {};
For[p = 1, p <= n, p++,
For[pi = 1, pi <= n, pi++,
If[RegionDimension[
RegionIntersection[MeshPrimitives[Last[vorl], 2][[p]],
Point[vlc[[1]][[pi]]]]] =!= -Infinity,
val1 = Append[val1, val[[1, pi]]]
]]];
For[i = 1, i <= T, i++,
sca = (i 0.01 + 1);
vl = Append[vl, rel[Last[vl], 1.05]];
vorl = Append[vorl, VorR[Last[vl]*sca, sca]];
vlc = Append[vlc,
Map[RegionCentroid, MeshPrimitives[Last[vorl], 2]]];
nsol = {};
dsol = {};
For[p = 1, p <= n, p++,
If[RegionDimension[
RegionIntersection[MeshPrimitives[Last[vorl], 2][[p]],
Point[vlc[[i]][[p]]]]] =!= -Infinity,
nsol = Append[nsol, val[[i, p]]],
For[pi = 1, pi <= n, pi++,
If[RegionDimension[
RegionIntersection[MeshPrimitives[Last[vorl], 2][[p]],
Point[vlc[[i]][[pi]]]]] =!= -Infinity,
nsol = Append[nsol, val[[i, pi]]]
]]]];
val = Append[val, nsol]
];
colr = ReplacePart[val, 1 -> val1];
cels = Table[
Table[Graphics[{RGBColor[0.5, 0.65, 0.5, colr[[j, p]]],
MeshPrimitives[vorl[[j]], 2][[p]]}], {p, n}], {j, T}];
ListAnimate[Table[Show[cels[[j]]], {j, T}]]
leading to
My previous code is even more clear if instead you take these values as
val = {Range[n]/n};
to get
This, however, is not ideal, because it is based on the previous cell centroid, and whether it intersects the current cell in the current mesh. In other words, all I'm doing is the following: I test whether the previous centroid intersects the corresponding new cell (in the list ordering). If yes, I don't change the order, and if not I track which centroid is contained in this cell and reorder the values of each cell colour (val) accordingly. However, if the cell movement is too quick (per time iteration), the previous centroids are not guaranteed to intersect the new cells and order is compromised. So this is one problem
I wonder, then, if there is both a neater and more clever way of doing this, I wish Mathematica had this in mind, because even when defining, for example, adjacency matrices over a system of ODEs that describe some dynamics between cells, if the mesh is not static, then ordering becomes an recurrent issue, and even adding or removing cells (vertices) in the mesh (adjacency graph) requires careful management in order to get the correct ordering. I also noticed that reordering seems to take place whenever the number of neighbours of cells change (or degree of the corresponding graph vertex).
Any ideas? Sorry for the long post.
Note: the code doesn't always work due to Voronoi seeds being away from the cropping region, but it doesn't interfere with the main point of the question. Simply run it again. This was written in Mathematica 12.1.
VoronoiMesh. $\endgroup$ – Alex Trounev Jun 30 at 0:28VoronoiMesh, otherwise there is no point in making this question. $\endgroup$ – sam wolfe Jun 30 at 1:05VoronoiMeshand Lloyd's relaxation algorithm. Do you try to improve this algorithm or just this code? $\endgroup$ – Alex Trounev Jun 30 at 12:37VoronoiMesh? $\endgroup$ – flinty Jun 30 at 12:58