I need two compiled helper functions:
getEdgesFromPolygons = Compile[{{f, _Integer, 1}},
Table[
{
Min[Compile`GetElement[f, i], Compile`GetElement[f, Mod[i + 1, Length[f], 1]]],
Max[Compile`GetElement[f, i], Compile`GetElement[f, Mod[i + 1, Length[f], 1]]]
},
{i, 1, Length[f]}
],
RuntimeAttributes -> {Listable},
Parallelization -> True,
CompilationTarget -> "C"
];
takeThread = Compile[{{data, _Integer, 1}, {ran, _Integer, 1}},
Part[data, ran[[1]] ;; ran[[2]]],
RuntimeAttributes -> {Listable},
Parallelization -> True,
CompilationTarget -> "C"
];
The following functions takes a `MeshRegion` and finds all pairs of neighboring two-dimensional `MeshCells`.
getPolygonPolygonAdjacencyList[R_MeshRegion] :=
Module[{pts, polygons, edgesfrompolygons, edges, edgelookupcontainer,
polyranges, signedpolygonsneighedges, polygonedgeadjacencymatrix},
pts = MeshCoordinates[R];
polygons = Developer`ToPackedArray[MeshCells[R, 2][[All, 1]]];
edgesfrompolygons =
Developer`ToPackedArray[Flatten[getEdgesFromPolygons[polygons], 1]];
edges = DeleteDuplicates[edgesfrompolygons];
edgelookupcontainer = SparseArray[
Rule[
Join[edges, Transpose[Transpose[edges][[{2, 1}]]]],
Join[Range[1, Length[edges]], -Range[1, Length[edges]]]
],
{Length[pts], Length[pts]}
];
polyranges = With[{acc = Accumulate[Length /@ polygons]},
Transpose[{Prepend[Most[acc] + 1, 1], acc}]];
signedpolygonsneighedges = takeThread[
Extract[
edgelookupcontainer,
edgesfrompolygons
],
polyranges
];
polygonedgeadjacencymatrix = With[{
n = Length[edges], m = Length[polygons],
data = Flatten[SortBy[#, Abs] & /@ signedpolygonsneighedges]
},
SparseArray @@ {Automatic, {m, n}, 0, {1, {
Developer`ToPackedArray[Join[{0}, Accumulate[Length /@ polygons]]],
Transpose[{Abs[data]}]
},
Sign[data]}}
];
Select[(Transpose[polygonedgeadjacencymatrix]["AdjacencyLists"]),
Length[#] == 2 &]
]
Testing with OP's example:
SeedRandom[123]
pts = RandomReal[1, {10, 2}];
R = VoronoiMesh[pts, MeshCellLabel -> {2 -> "Index"}]
Graph[
UndirectedEdge @@@ getPolygonPolygonAdjacencyList[R],
VertexLabels -> "Name"
]
[![enter image description here][1]][1]
[![enter image description here][2]][2]
Speed test
SeedRandom[123]
pts = RandomReal[1, {10000, 2}];
R = VoronoiMesh[pts, MeshCellLabel -> {2 -> "Index"}]; // AbsoluteTiming
getPolygonPolygonAdjacencyList[R]; // AbsoluteTiming
> {0.639544, Null}
>
> {0.137034, Null}
[1]: https://i.sstatic.net/e1QhF.png
[2]: https://i.sstatic.net/7p6XR.png