Does this do what you have in mind?
pts = {{1, 1}, {1, 2}, {1234512345, 1234512345}, {1234512345, 1234512346}};
edges = Join @@ MapIndexed[
Thread[{#1, #2[[1]]}] &,
Nearest[
pts -> "Index",
pts,
{\[Infinity], 1},
DistanceFunction -> ManhattanDistance
][[All, 2 ;;]]
];
SparseArray`StronglyConnectedComponents[SparseArray[edges ->_]]
For diagonal connectivity the same code with ManhattanDistance
replaced by ChessboardDistance
should work.
Edit
At least for the nondiagonal connectivity, this can be sped up by a factor of ten by avoiding Nearest
. We can do so because we can find out the horizontal neighbors by going through the rows of sparse image matrix. This is easy to do for a matrix in CSR (compressed sparse row) format. For finding out the vertical neighbors, we just do that for the transposed image, too.
First a couple of helper functions:
QuickSparseArray[rp_?VectorQ, ci_?MatrixQ, vals_?VectorQ, dims_?VectorQ, background_ : 0] :=
With[{data = {Automatic, dims, background, {1, {rp, ci}, vals}}},
SparseArray @@ data
];
ThreadCount[] := "ParallelThreadNumber" /. ("ParallelOptions" /. SystemOptions["ParallelOptions"]);
JobPointers[jobCount_Integer?Positive, threadCount_Integer?Positive] :=
Ceiling[Subdivide[0, jobCount, Min[threadCount, jobCount]]];
cFindHorizontalNeighbors =
Compile[{{rp, _Integer, 1}, {ci, _Integer, 2}, {idx, _Integer,
1}, {start, _Integer}, {end, _Integer}},
Block[{bag, col, nextcol, i, j},
bag = Internal`Bag[Most[{0}]];
Do[
Do[
col = Compile`GetElement[ci, k, 1];
nextcol = Compile`GetElement[ci, k + 1, 1];
If[nextcol == col + 1,
i = Compile`GetElement[idx, k];
j = Compile`GetElement[idx, k + 1];
Internal`StuffBag[bag, i];
Internal`StuffBag[bag, j];
Internal`StuffBag[bag, j];
Internal`StuffBag[bag, i];
]
, {k, Compile`GetElement[rp, row] + 1,
Compile`GetElement[rp, row + 1] - 1}]
, {row, start, end}]; Partition[Internal`BagPart[bag, All], 2]
],
CompilationTarget -> "C",
RuntimeAttributes -> {Listable},
Parallelization -> True,
RuntimeOptions -> "Speed"
];
Now the actual computations:
n = Length[A["NonzeroValues"]];
B = QuickSparseArray[A["RowPointers"], A["ColumnIndices"], Range[n], Dimensions[A]];
BT = Transpose[B];
Bptr = JobPointers[Length[B], ThreadCount[]];
BTptr = JobPointers[Length[BT], ThreadCount[]];
edges = Join[
Join @@ cFindHorizontalNeighbors[B["RowPointers"], B["ColumnIndices"], B["NonzeroValues"], Most[Bptr] + 1, Rest[Bptr]],
Join @@ cFindHorizontalNeighbors[BT["RowPointers"], BT["ColumnIndices"], BT["NonzeroValues"], Most[BTptr] + 1, Rest[BTptr]]
];
components = SparseArray`StronglyConnectedComponents[SparseArray[edges -> _, {n, n}]];
Bonus: We can get the colored sparse array with
colors = Normal[SparseArray[ Join @@ components -> Join @@ (Range[Length[components]] Unitize[components]), n]];
Acolored = QuickSparseArray[A["RowPointers"], A["ColumnIndices"], colors, Dimensions[A]];
Now we can execute Colorize[Acolored]
and behold.
In principle, one could do a similar thing for the diagonal lookup. I am just not in the modd for that...