# Finding the size of connected clusters with specific patterns in an array

I have a $w \times h$ array with either $-1$ or $+1$.

I wish to monitor the size of the largest connected clusters with specific patterns that occur during a Monte Carlo simulation, i.e.,

• a collection of $+1$'s,
• a checkerboard pattern of $-1$ and $+1$'s, and,
• horizontal and vertical stripes of $-1$ and $+1$'s.

Moreover, periodic boundary conditions apply such that opposite boundaries are connected.

Example arrays (checkerboard, stripes) can be imported by,

Uncompress[Import["https://pastebin.com/raw/0GaUtiFy"]]
Uncompress[Import["https://pastebin.com/raw/uF5UagqU"]]


Ideally I am looking for a fast solution, such that monitoring every iteration of the Monte Carlo simulation will not slow it down to a crawl.

The only idea I have had so far was to use ComponentMeasurements but was not able to specify the patterns. Any help would be greatly appreciated.

## 1 Answer

This might not be the fasted method but it will get you started:

Preparation: The working horse is getComponent while periodicBoundaryCorrection only accounts for the fact that MorphologicalComponents does not allow for periodic boundary conditions. See the also here for a more detailed description.

checkermask = Array[Mod[Plus[##], 2] &, {3, 3}];
hormask = Array[Mod[#1, 2] &, {3, 3}];
vermask = Array[Mod[#2, 2] &, {3, 3}];

ClearAll[periodicBoundaryCorrection];
periodicBoundaryCorrection[
A_?MatrixQ,
OptionsPattern[{CornerNeighbors -> True}]
] := Module[{maxA, a, b, c, d, pos, r1, r2, r3, r4, r5, r6, edges, α, β, colorcomp, cols},
maxA = Max[A];

a = A[[1]];
b = A[[-1]];
pos = DeleteCases[Range[Length[a]] Unitize[a b], 0];
r1 = Sort /@ Transpose[{a[[pos]], b[[pos]]}];

c = A[[All, 1]];
d = A[[All, -1]];
pos = DeleteCases[Range[Length[c]] Unitize[c d], 0];
r2 = Sort /@ Transpose[{c[[pos]], d[[pos]]}];

edges = Union[r1, r2];

If[OptionValue[CornerNeighbors],

pos = DeleteCases[Range[Length[a] - 1] Unitize[Rest[a] Most[b]], 0];
r1 = Sort /@ Transpose[{a[[pos + 1]], b[[pos]]}];
pos = DeleteCases[Range[Length[b] - 1] Unitize[Rest[b] Most[a]], 0];
r2 = Sort /@ Transpose[{b[[pos + 1]], a[[pos]]}];

pos = DeleteCases[Range[Length[c] - 1] Unitize[Rest[c] Most[d]], 0];
r3 = Sort /@ Transpose[{c[[pos + 1]], d[[pos]]}];
pos = DeleteCases[Range[Length[d] - 1] Unitize[Rest[d] Most[c]], 0];
r4 = Sort /@ Transpose[{d[[pos + 1]], c[[pos]]}];

α = c[[1]];
β = d[[-1]];
r5 = If[α β != 0, {Sort[{α, β}]}, {}];
α = c[[-1]];
β = d[[1]];
r6 = If[α β != 0, {Sort[{α, β}]}, {}];

edges = Union[edges, r1, r2, r3, r4, r5, r6]
];
edges++;
If[Length[edges] == 0,
A,
colorcomp = SparseArrayStronglyConnectedComponents[
SparseArray[
Join[edges, Transpose[Transpose[edges][[{2, 1}]]]] ->
1, {maxA + 1, maxA + 1}, 0]
];
cols = Compile[{{idx, _Integer, 1}, {acc, _Integer, 1}},
Block[{colors, j, threshold},
colors = Table[0, {i, 1, acc[[-1]]}];
j = 0;
threshold = CompileGetElement[acc, j + 1];
colors[[idx]] = Table[
If[i > threshold,
j++;
threshold = CompileGetElement[acc, j + 1];
];
j, {i, 1, Length[idx]}];
colors
]
][
Join @@ colorcomp,
Accumulate[Length /@ colorcomp]
];
Compile[{{a, _Integer, 1}, {cols, _Integer, 1}},
cols[[a + 1]],
RuntimeAttributes -> {Listable},
Parallelization -> True
][A, cols]
]
];

getComponent[A_?MatrixQ, mask_] := periodicBoundaryCorrection[
MorphologicalComponents[
With[{img = Binarize[Image[A], 0]},
With[{
a = ImageCorrelate[img, mask/Total[mask, 2], Padding -> "Periodic"],
b = ImageCorrelate[img, (1 - mask)/Total[1 - mask, 2],
Padding -> "Periodic"]
},
With[{ab = ImageMultiply[a, b]},
Binarize[ImageAdd[ImageSubtract[a, ab], ImageSubtract[b, ab]],
0.99]
]
]
]
]
]


Application:

A = Uncompress[Import["https://pastebin.com/raw/0GaUtiFy"]];
B = Uncompress[Import["https://pastebin.com/raw/uF5UagqU"]];

GraphicsGrid[{
{
Image[A],
Colorize[getComponent[A, checkermask]],
Colorize[getComponent[A, hormask]],
Colorize[getComponent[A, vermask]]
},
{
Image[B],
Colorize[getComponent[B, checkermask]],
Colorize[getComponent[B, hormask]],
Colorize[getComponent[B, vermask]]
}
}, ImageSize -> Large]


• Thank you, this definitely got me started. Reverse engineering your solution I found that getComponent[A_, mask_] := With[{B = -A ListCorrelate[mask, A, {2, 2}]}, ImageData[Binarize[Image[B], 3]]] almost does the trick. The only thing left is to implement the boundary conditions Apr 11, 2018 at 19:30
• Ah sorry, I missed that. You have to set Padding - > "Periodic" in ListCorrelate. At least that should make the input for MorphologicalComponents correct. MorphologicalComponents won't match opposing components correctly. It has a Padding options but that does not allow for "Periodic"`... Apr 11, 2018 at 19:32
• No sorry I will clarify, the components should found with periodic boundary conditions, i.e. in the checkerboard example there should only be one cluster found, instead of four Apr 11, 2018 at 19:40
• In any case this question already has a answer Apr 11, 2018 at 19:52
• Thanks for the hint. It was too late; already added a fix. Apr 11, 2018 at 20:19