Symmetry-Reduced Enumeration and Canonical Forms for Grid-Based Density Clustering Under the Hyperoctahedral Group

1. Introduction

Clustering is the unsupervised partition of a dataset into groups of similar points; density-based clustering defines those groups as connected regions of high density separated by sparse ones. DBSCAN [1] and its ordering-based refinement OPTICS [2] are the canonical representatives, widely used because they find clusters of arbitrary shape without fixing their number in advance. Their principal cost is the neighbourhood query, and a standard remedy is to impose a regular grid on the feature space, aggregate points into axis-aligned bins, and run the density logic on bins rather than on individual points. FlowGrid [3] and FLOPTICS [4,5,6] are recent members of this grid-based density family; FLOPTICS, in particular, was introduced for the automated gating of flow-cytometry data—a setting in which panels typically measure on the order of ten to forty channels, though gating is usually performed on low-dimensional projections of two to four channels at a time, precisely the regime where $\left|B_d\right|\in \{8,48,384\}$ is small. Axis relabelling is not a hypothetical symmetry here: the order in which channels are written to a data file, and the orientation conventions of a given instrument or compensation step, are recording choices rather than biological facts, so two laboratories analysing the same sample may present grids differing exactly by an element of $B_d$. FLOPTICS clusters bins produced by an equal-width partition of each measured channel, and belongs to a family of automated gating methods that includes FLOCK [7], flowPeaks [8], and FlowSOM [9].

A grid is a highly symmetric object. Permuting the $d$ feature axes and independently reflecting any axis sends bins to bins and leaves every notion of grid adjacency intact. The collection of these operations is a finite group—the hyperoctahedral group $B_d=C_2\wr S_d$ of order $2^{\hspace{0.17em}d}d!$, classically the symmetry group of the $d$-cube and of the integer lattice ${\mathbb{Z}}^d$, and the natural symmetry acting on functions over a product space of the kind studied in ${\left[m\right]}^d$ [10]; the underlying integer-lattice geometry is standard [11]. Although grid-based clustering algorithms operate directly on this symmetric object, they treat the axes as ordered and oriented, so the symmetry is present in the problem but absent from its analysis. The purpose of this paper is to make the symmetry explicit and to derive its exact consequences for the enumeration, canonicalisation, and reproducibility of grid-clustering instances.

Our motivation is concrete and recent. In a companion study of the DR Code—a two-dimensional barcode whose $3\times 3$ block layout admits XOR-based strip recovery [12]—the authors showed [13] that the recoverable layouts form the orbit space of an order-72 group and that passing to orbit representatives compresses the layout catalogue and yields a canonical-form deduplication scheme. The geometric part of that group is exactly $D_4$, the symmetry group of the square. A companion study extends the code to a $3\times 3\times 3$ logical cube, where the corresponding symmetry group is $\left(C_3\times C_3\times C_3\right)⋊O_h$ of order $1,296$ and the geometric factor is the full octahedral group $O_h$ of the cube [14]. Both are instances of the hyperoctahedral family: $D_4=B_2$ and $O_h=B_3$. The present paper takes the natural next step by removing the barcode-specific recoverability constraint and studying the pure grid symmetry $B_d$ in the setting where the grid hosts a clustering instance rather than a constrained code. In doing so we move from a finite constrained catalogue to the unconstrained—and unbounded—space of bin-occupancy functions, which changes both the enumeration problem and the meaning of a canonical form.

We are careful about what symmetry buys. A generic occupancy function has trivial stabiliser, so its orbit has full size and no bins are forced to coincide; we therefore make no per-run speedup claim. The orbit structure is instead useful across instances—canonical indexing, deduplication of symmetry-closed libraries, and caching across analyses that present the same measurement under a different axis convention—and we quantify these benefits.

A second motivation is reproducibility. Grid-based clustering pipelines are increasingly deployed as standardised tools—FLOPTICS, for instance, is intended as a reproducible alternative to subjective manual gating—and reproducibility requires that the output not depend on incidental choices such as the order in which feature channels happen to be listed or the sign convention of an axis. Equivariance under $B_d$ is the exact statement that these choices are incidental: two operators who list the channels in different orders, or orient an axis oppositely, obtain clusterings differing by precisely the corresponding relabelling. A canonical axis order removes this degree of freedom and makes the pipeline deterministic with respect to axis presentation. This is a modest but genuine guarantee, and it falls out of the same group action that governs the enumeration and caching results.

Two questions organise what follows. (Q1) In what precise sense is grid-based density clustering independent of how the feature axes happen to be ordered and oriented, and what is the group that expresses this? (Q2) Given that group, how many genuinely distinct clustering instances are there on a grid, and can each be given a canonical name cheaply enough to be useful? Q1 is answered by the equivariance theorem, Q2 by the enumeration and canonical-form results.

The contributions of this paper are the following. (C1) We formalise the action of $B_d$ on bin-occupancy functions and prove that it preserves the Chebyshev bin-adjacency that underlies grid density reachability (Section 4, Theorem 1). (C2) We prove that grid-based density clustering is $B_d$-equivariant: relabelling the axes relabels the clusters and nothing more (Section 5, Theorem 2). (C3) We give a Burnside/cycle-index count of bin-occupancy configurations up to $B_d$, specialise the classical signed-cycle-type computation of the cycle index of $B_d$ to the grid ${\left[m\right]}^d$; the resulting closed form for the octahedral case, valid at every resolution and split by the parity of the resolution, is our main enumerative contribution and we have not found it recorded elsewhere. The dihedral case ($d=2$) is derived the same way as a worked illustration, and we verify all of it against exhaustive orbit enumeration (Section 6, Theorems 3 and 4, Corollary 2). (C4) We present a canonical-form algorithm (Algorithm 1), prove its correctness and $O\left(2^{\hspace{0.17em}d}d!\hspace{0.17em}{m}^{\hspace{0.17em}d}\right)$ complexity, and delimit its honest range of usefulness (Section 7, Theorem 5). Throughout, two supporting propositions record the orbit-size structure (Proposition 1) and the asymptotic orbit count (Proposition 2), and a fully worked $2\times 2$ micro-example makes every ingredient checkable by hand. Section 2 reviews related work, Section 3 fixes notation, Section 8 discusses scope and applications, and Section 9 concludes. In total we prove five theorems, one lemma, two propositions and four corollaries, and give two algorithms. Every quantitative claim in the paper was checked computationally using the open-source scripts released with the paper—exhaustively on the relevant finite structures wherever they are small enough to enumerate, and by large randomised trials elsewhere.

2. Background and Related Work

2.1. Density-Based and Grid-Based Clustering

DBSCAN defines clusters through core points, whose $\epsilon$-neighbourhood contains at least a threshold number of points, and the density-reachability relation they induce. OPTICS replaces the single global $\epsilon$ by an ordering of the points annotated with reachability distances, from which clusterings at multiple density levels can be extracted [2,15]. Both incur neighbourhood queries whose naïve cost is quadratic and whose spatially indexed cost is $O\left(n\mathrm{l}\mathrm{o}\mathrm{g}n\right)$ [16]. Grid-based methods sidestep the query by bucketing points into a regular array of bins: FlowGrid and FLOPTICS partition each dimension into equal-width intervals, keep only non-empty bins, and apply density logic to bin counts and bin adjacency; this generalised density-based view, in which clusters are connected dense cells of a space partition, goes back to GDBSCAN [17]. FLOPTICS specifically runs an OPTICS ordering over the occupied bins, so the clusters it extracts are determined by bin counts and the Chebyshev adjacency of the grid—the two ingredients our model retains—while discarding sub-bin coordinates [4,5,6]. This reduces the working set from $n$ points to the number of occupied bins and replaces metric neighbourhood tests by index arithmetic. The methods we analyse are exactly of this grid-based density type [3,4,5,18]; our results apply to any algorithm whose clusters are determined by bin occupancy and grid adjacency.

Table 1 lists the methods in this family and what our model retains from each.

To make the object of study concrete, Figure 1 illustrates the grid-based density-clustering pipeline on a three-cluster synthetic example in the style of automated flow-cytometry gating. A reachability ordering (panel a) exposes the cluster structure as valleys; the feature space is discretised into an axis-aligned grid whose non-empty bins are the working set (panel b); and the clusters recovered on that grid correspond to the dense regions of the original data (panel c). The figure is schematic and carries no experimental claim—the empirical evaluation of FLOPTICS itself is reported in its source work [4,5,6]—but it fixes the picture our analysis depends on: everything downstream is determined by which bins are occupied and how they are adjacent, and it is precisely this grid, panel (b), on which the hyperoctahedral group $B_d$ acts in the remainder of the paper.

2.2. Symmetry of the Grid and the Hyperoctahedral Group

The symmetry group of the hypercubic lattice ${\mathbb{Z}}^d$—equivalently of the $d$-cube—is the hyperoctahedral group of all signed permutation matrices, $B_d=C_2^{\hspace{0.17em}d}⋊S_d$, of order $2^{\hspace{0.17em}d}d!$. Its elements are the signed permutations: a permutation of the $d$ axes composed with an independent choice of reflection on each axis. In the combinatorial literature $B_d$ is the group of the “signed” or “hyperoctahedral” symmetric functions, and its cycle index is the classical tool for enumerating structures on the cube and lattice [19,20,21]. We use $B_d$ in precisely this role, but applied to the array of clustering bins rather than to the cube’s vertices.

2.3. Orbit Enumeration by Burnside’s Lemma

Counting configurations up to a group action is governed by the Cauchy–Frobenius–Burnside lemma: the number of orbits equals the average number of configurations fixed by a group element [22,23]. When configurations are colourings of a set acted on by a group, the count is delivered by Pólya’s cycle-index theorem [24,25]; general treatments appear in [21,26,27]. These techniques have been applied to lattice-square designs [28], whose isomorphism group is built from a hyperoctahedral and a dihedral factor and whose classes are counted by Burnside; to the enumeration of polyomino symmetry classes on the square lattice under $C_4$ and $D_4$ [29]; and to the classification of full-dimensional $0/1$-polytopes through the cycle index of $B_n$ [30]. Our enumeration (Theorem 3) is a Burnside average over $B_d$ acting on occupancy functions; it is methodologically classical, and our contribution is its transfer to, and interpretation for, clustering instances rather than any new counting machinery.

2.4. Canonical Forms and Reproducibility

A canonical form assigns to each object a unique representative of its equivalence class, so that two objects are equivalent exactly when their canonical forms coincide. Canonical labelling of graphs (the individualisation–refinement algorithms behind nauty and Traces) [31] and canonical string encodings of molecules (canonical SMILES) [32] are the standard instances; both underpin database deduplication, indexing, and reproducible pipelines. Our canonical form (Theorem 5) is of the same spirit but for a far smaller and fully explicit group: because $\left|B_d\right|=2^{\hspace{0.17em}d}d!$ is fixed and independent of the data, canonicalisation is a bounded loop over group elements rather than a search over an exponential automorphism space, and its correctness is immediate. This is the practical contrast with general graph canonicalisation: the individualisation–refinement search behind nauty and Traces may, in adversarial cases, backtrack through exponentially many candidate labellings, whereas the canonical form used here performs exactly $2^{\hspace{0.17em}d}d!$ push-forwards with no backtracking and no data-dependent branching, so its cost is a fixed constant ($8$, $48$, $384$ for $d=2,3,4$) times the grid size.

2.5. Relation to the DR Code

The recoverable layouts of the $3\times 3$ DR Code form the orbit space of a group $\left(C_3\times C_3\right)⋊D_4$ of order 72, whose geometric factor $D_4$ acts by rotating and reflecting the square grid while the cyclic factor permutes row and column labels subject to the recoverability constraint [13]. The present paper isolates the geometric factor and generalises it: we drop the recoverability constraint, replace the fixed $3\times 3$ board by an arbitrary $m^{\hspace{0.17em}d}$ grid, and replace $D_4=B_2$ by the full hyperoctahedral group $B_d$. Two differences are worth stating plainly. First, the DR Code group mixes a geometric ($D_4$) and a combinatorial ($C_3\times C_3$) factor tied to its recoverability rule; our group is purely the geometric grid symmetry $B_d$, appropriate because a clustering grid carries no analogue of the row/column-label constraint. Second, the DR Code acts on a finite constrained set of 2,592 layouts, whereas $B_d$ acts on the unbounded space of occupancy functions, so our enumeration is parameterised by grid resolution and occupancy alphabet rather than being a single finite number. The DR Code is thus the $d=2$, constrained special case that motivates the unconstrained, general-$d$ theory developed here.

Table 2 situates the present work among prior applications of hyperoctahedral and lattice symmetry. The machinery—signed permutations with Burnside/cycle-index counting—is shared; what is new is the object and its purpose: an unbounded space of clustering instances analysed for canonical indexing and equivariant reproducibility rather than a finite classification.

2.6. Relation to Equivariant Machine Learning

The idea that a model should respect the symmetries of its input is central in modern machine learning: group-equivariant convolutional networks [35] build layers that commute with a symmetry group, and the geometric-deep-learning programme [36] organises architectures around this principle, with grids as the prototypical example. Theorem 2 is a statement of the same kind for a classical clustering procedure rather than a learned network, with two differences: our grids are finite and unpadded, so the relevant group is the finite point group $B_d$ rather than a translation group; and the payoff here is enumerative rather than sample efficiency. A third connection is data augmentation—a $B_d$-closed library is exactly an axis-augmented suite, so the deduplication of Section 8 says such a suite carries no more information than its canonical reduction. Hierarchical methods such as HDBSCAN [15] share the same structure on a binned space, so the results apply to them unchanged.

3. Preliminaries

For convenience, Table 3 collects the notation used throughout.

Fix a dimension $d\ge 1$ and a per-axis resolution $m\ge 2$. The grid is the set of bins

$$\Gamma \hspace{0.33em}=\hspace{0.33em}{\left[m\right]}^d\hspace{0.33em}=\hspace{0.33em}\{0,1,\dots ,m-1{\}}^d, \left|\Gamma \right|=m^{\hspace{0.17em}d}.$$

Two bins $u,v\in \Gamma$ are adjacent, written $u\sim v$, if they differ and are within Chebyshev distance one, i.e. $u\ne v$ and ${\mathrm{m}\mathrm{a}\mathrm{x}}_k\left|u_k-v_k\right|\le 1$. This is the king-move adjacency used by grid-based density clustering to decide when two occupied bins belong to the same dense region.

An occupancy function is a map $w:\Gamma \to A$ into a finite alphabet $A$ of size $q=\left|A\right|$; $w\left(b\right)$ records the point count of bin $b$, optionally quantised into $q$ density levels, with $0$ denoting an empty bin. The set of occupancy functions is $A^{\Gamma }$, of size $q^{\hspace{0.17em}{m}^d}$.

It is convenient to record two special cases immediately. At $d=2$ the group $B_2=C_2^2⋊S_2$ has order $8$ and is the dihedral group $D_4$ of the square: the two axis reflections and the axis swap generate the four rotations and four reflections of a square grid. At $d=3$ the group $B_3=C_2^3⋊S_3$ has order $48$ and is the full octahedral group $O_h$ of the cube. These are exactly the geometric symmetry groups appearing in the 2D and 3D DR Code, which is the sense in which the hyperoctahedral action generalises that setting to arbitrary dimension.

We also fix the cycle-index notation used in Section 6. For a permutation $\sigma$ of a finite set $X$ write its cycle type as the formal monomial $\prod _{k\ge 1}{a}_k^{\hspace{0.17em}{j}_k\left(\sigma \right)}$, where $j_k\left(\sigma \right)$ is the number of $k$-cycles of $\sigma$. The cycle index of a group $H$ acting on $X$ is $Z_H^X\left(a_1,a_2,\dots \right)={\displaystyle \frac{1}{\left|H\right|}}\sum _{\sigma \in H}\prod _k{a}_k^{\hspace{0.17em}{j}_k\left(\sigma \right)}$. Setting every $a_k=q$ evaluates the cycle index to ${\displaystyle \frac{1}{\left|H\right|}}\sum _{\sigma }{q}^{\hspace{0.17em}c\left(\sigma \right)}$ with $c\left(\sigma \right)=\sum _k{j}_k\left(\sigma \right)$ the total number of cycles; by Theorem 3 this is the number of $q$-colourings of $X$ up to $H$.

The hyperoctahedral group is

$$B_d\hspace{0.33em}=\hspace{0.33em}{C}_2\wr S_d\hspace{0.33em}=\hspace{0.33em}{C}_2^{\hspace{0.17em}d}⋊S_d, \left|B_d\right|=2^{\hspace{0.17em}d}\hspace{0.17em}d!,$$

realised as signed permutations. We write an element as a pair $g=\left(\pi ,\epsilon \right)$ with $\pi \in S_d$ a permutation of axes and $\epsilon \in \{+1,-1{\}}^d$ a per-axis reflection sign. Composition follows the wreath-product rule: for $g_1=\left({\pi }_1,{\epsilon }_1\right)$ and $g_2=\left({\pi }_2,{\epsilon }_2\right)$ one has $g_1{g}_2=\left({\pi }_1{\pi }_2,\epsilon \right)$ with ${\epsilon }_i={\epsilon }_{2,i}\hspace{0.17em}{\epsilon }_{1,{\pi }_2\left(i\right)}$, and the identity is $\left(\mathrm{i}\mathrm{d},\left(+1,\dots ,+1\right)\right)$. We use two standard facts about a finite group $H$ acting on a set: the orbit of a point $x$ is $H\cdot x=\{h\cdot x:h\in H\}$, its stabiliser is $H_x=\{h:h\cdot x=x\}$, and the orbit–stabiliser theorem gives $\left|H\cdot x\right|=\left|H\right|/\left|H_x\right|$; and the number of orbits is the Burnside average of fixed-point counts. Both are applied below to $H=B_d$ acting on occupancy functions.

4. The Group Action on Grid Bins

We first make $B_d$ act on bins, then on occupancy functions, and verify the action respects grid adjacency—the property that makes it compatible with density reachability. The three generating operations are shown schematically in Figure 2.

Definition 1 (action on bins). For $g=\left(\pi ,\epsilon \right)\in B_d$ and a bin $x=\left(x_0,\dots ,x_{d-1}\right)\in \Gamma$ , define $g\cdot x=y$ by

$$y_{\pi \left(i\right)}\hspace{0.33em}=\hspace{0.33em}\left\{\begin{array}{ll}{x}_i,& {\epsilon }_i=+1,\\ \left(m-1\right)-x_i,& {\epsilon }_i=-1,\end{array}\right. i=0,\dots ,d-1.$$

Thus $g$ permutes the axes by $\pi$ and reflects the selected axes. The induced action on occupancy functions is by push-forward, $\left(g\cdot w\right)\left(b\right)=w\left(g^{-1}\cdot b\right)$.

Theorem 1 (well-defined adjacency-preserving action). The map of Definition 1 is a left action of $B_d$ on $\Gamma$ by bijections, and it preserves adjacency: $u\sim v\iff g\cdot u\sim g\cdot v$ for every $g\in B_d$ . Consequently the push-forward is a left action of $B_d$ on $A^{\Gamma }$ , and it preserves the bin-adjacency graph on which grid density reachability is defined. 

Proof. Each single-axis reflection $x_i\mapsto \left(m-1\right)-x_i$ is an involution of $\{0,\dots ,m-1\}$, and axis permutation by $\pi$ is a bijection of coordinate positions; their combination is therefore a bijection of $\Gamma$, with inverse ${\left(\pi ,\epsilon \right)}^{-1}=\left({\pi }^{-1},\epsilon \mathrm{\text{'}}\right)$ where $\epsilon {\mathrm{\text{'}}}_{\pi \left(i\right)}={\epsilon }_i$. The identity $\left(\mathrm{i}\mathrm{d},\left(+1,\dots ,+1\right)\right)$ acts as the identity map, and a direct computation shows $\left(g_1{g}_2\right)\cdot x=g_1\cdot \left(g_2\cdot x\right)$, so the map is a left action. For adjacency, write ${\delta }_k=\left|u_k-v_k\right|$. Axis permutation only relabels the coordinates, so it permutes the multiset $\{{\delta }_0,\dots ,{\delta }_{d-1}\}$ and leaves its maximum unchanged. A reflection on axis $i$ replaces $u_i,v_i$ by $\left(m-1\right)-u_i,\left(m-1\right)-v_i$, whose absolute difference is again ${\delta }_i$. Hence ${\mathrm{m}\mathrm{a}\mathrm{x}}_k{\delta }_k$ is invariant and $u\ne v\iff g\cdot u\ne g\cdot v$, so $u\sim v\iff g\cdot u\sim g\cdot v$. The push-forward inherits the action axioms from the action on bins, and because adjacency is preserved, the bin-adjacency graph is mapped isomorphically to itself. $◼$

We verified Theorem 1 exhaustively: for $d\in \{2,3,4\}$ and several resolutions $m$, all $\left|B_d\right|$ elements act as bijections of $\Gamma$, the wreath-product composition rule agrees with the composition of the corresponding bin-maps on every bin, and every adjacent pair maps to an adjacent pair. At $d=2$ the group has order $8$ and coincides with $D_4$; at $d=3$ it has order $48$ and coincides with $O_h$—precisely the geometric groups of the 2D and 3D DR Code (Table 4).

To make the action concrete, take $d=2$ and $m=3$, and index the nine bins $0,\dots ,8$ in row-major order so that bin $\left(i,j\right)$ receives index $3i+j$. The axis swap $\pi =\left(0\hspace{0.33em}1\right)$ with trivial signs sends $\left(i,j\right)\mapsto \left(j,i\right)$; as a permutation of bins it is the transposition pattern fixing the diagonal $\{0,4,8\}$ and swapping $1\leftrightarrow 3$, $2\leftrightarrow 6$, $5\leftrightarrow 7$, i.e. a matrix transpose. The reflection of axis $0$ with trivial permutation sends $\left(i,j\right)\mapsto \left(2-i,j\right)$, exchanging the top and bottom rows. Composing the eight signed permutations of $B_2$ reproduces precisely the eight symmetries of the square acting on the $3\times 3$ board, confirming $B_2\cong D_4$ at the level of the bin permutation representation. The same construction at $d=3$, $m=3$ enumerates the $48$ bin permutations of $O_h$ on the $27$-cell cube.

5. Equivariance of Grid Density Clustering

We now show that clustering commutes with the group action, after fixing the model precisely.

Definition 2 (grid density clustering).Fix a density threshold$\tau \in A$. A bin$b\in \Gamma$is dense for the occupancy function$w$if$w\left(b\right)\ge \tau$; write$D_{\tau }\left(w\right)=\{b\in \Gamma :w\left(b\right)\ge \tau \}$for the set of dense bins. The grid density clustering of$w$is

$$\mathcal{C}\left(w\right)\hspace{0.33em}=\hspace{0.33em}\left\{\mathrm{connected}\mathrm{components}\mathrm{of}\mathrm{the}\mathrm{graph}\left(D_{\tau }\left(w\right),\sim \right)\right\},$$

a partition of $D_{\tau }\left(w\right)$ into clusters, where $\sim$ is the Chebyshev adjacency ofSection 3. This is the grid-based, single-threshold specialisation of density clustering that FlowGrid applies directly and that an OPTICS-based method such as FLOPTICS yields when clusters are extracted at a fixed reachability level; it is determined by the pair (which bins are dense, which dense bins are adjacent) and by nothing else. The intermediate reachability ordering produced by OPTICS is not itself claimed to be canonical: it may depend on processing order and tie-breaking, whereas the extracted clustering treated here does not. 

With Definition 2 in place, the equivariance statement is immediate to formulate.

Theorem 2 (equivariance).For every$g\in B_d$and every occupancy function$w$,

$$\mathcal{C}\left(g\cdot w\right)\hspace{0.33em}=\hspace{0.33em}g\cdot \mathcal{C}\left(w\right),$$

where$g$acts on a partition by acting on each of its bins. In particular the number of clusters and the multiset of cluster sizes are$B_d$-invariant, and the cluster containing bin$b$is mapped by$g$to the cluster containing$g\cdot b$. 

Proof. A bin $b$ is dense for $g\cdot w$ iff $\left(g\cdot w\right)\left(b\right)=w\left(g^{-1}\cdot b\right)\ge \tau$, i.e. iff $g^{-1}\cdot b$ is dense for $w$; thus the dense set of $g\cdot w$ is the $g$-image of the dense set of $w$. By Theorem 1 the map $b\mapsto g\cdot b$ is an adjacency graph isomorphism, and a graph isomorphism carries connected components to connected components bijectively. Therefore the components of the dense set of $g\cdot w$ are exactly the $g$-images of the components of the dense set of $w$, which is the claimed identity $\mathcal{C}\left(g\cdot w\right)=g\cdot \mathcal{C}\left(w\right)$. The stated invariants follow because a bijection preserves cardinalities. $◼$

Corollary 1 (invariance of clustering summaries).Every function of a clustering that depends only on the multiset of cluster sizes—the number of clusters, the largest cluster size, the number of clusters exceeding a given size, and any index or entropy computed from the size distribution—is a$B_d$ -invariant of the occupancy function, hence constant on orbits and tabulable per orbit rather than per instance. 

Proof. By Theorem 2 the clustering of $g\cdot w$ is the image of the clustering of $w$ under the bijection $b\mapsto g\cdot b$. A bijection preserves the multiset of block sizes, so any function of that multiset agrees on $w$ and $g\cdot w$ for every $g\in B_d$, i.e. is constant on the orbit. $◼$

Equivariance is the precise sense in which the ordering and orientation of the feature axes are irrelevant to the clustering: relabelling the axes relabels the output and changes nothing else. This has an immediate practical reading. A pipeline that fixes a canonical axis order (Section 7) produces identical clusters regardless of the order in which channels were originally recorded, removing a spurious source of run-to-run variation. We verified Theorem 2 on $6,000$ (occupancy, group-element, threshold) triples—$d=2$ with $m\in \{3,4,5\}$ and $d=3$ with $m\in \{3,4\}$, using randomised thresholds: in every case the partition of $g\cdot w$ equalled the $g$-image of the partition of $w$. Figure 3 depicts the commuting square that the theorem asserts.

6. Enumeration of Bin-Occupancy Classes

Because clustering is $B_d$-equivariant, two occupancy functions in the same orbit are indistinguishable up to axis relabelling and yield the same clustering up to that relabelling. It is therefore natural to count occupancy functions up to $B_d$. Let $N\left(d,m,q\right)$ denote the number of orbits of $B_d$ on $A^{\Gamma }$ with $\left|A\right|=q$.

Theorem 3 (Burnside enumeration).The number of bin-occupancy configurations on the$m^{\hspace{0.17em}d}$grid up to$B_d$is

$$N\left(d,m,q\right)\hspace{0.33em}=\hspace{0.33em}{\displaystyle \frac{1}{\left|B_d\right|}}\sum _{g\in B_d}{q}^{\hspace{0.17em}c\left(g\right)}\hspace{0.33em}=\hspace{0.33em}{\displaystyle \frac{1}{2^{\hspace{0.17em}d}d!}}\sum _{g\in B_d}{q}^{\hspace{0.17em}c\left(g\right)},$$

where$c\left(g\right)$is the number of cycles of$g$acting as a permutation of the$m^{\hspace{0.17em}d}$bins. Equivalently,$N\left(d,m,q\right)=Z_{B_d}^{\Gamma }\left(q,q,\dots ,q\right)$, the cycle index of the bin-permutation representation of$B_d$evaluated with every variable set to$q$. 

Proof. The group $B_d$ acts on the finite set $A^{\Gamma }$ of colourings of $\Gamma$ by $q$ colours through the push-forward. By the Cauchy–Frobenius–Burnside lemma the number of orbits equals the average number of fixed colourings, $N={\displaystyle \frac{1}{\left|B_d\right|}}\sum _g\left|\mathrm{F}\mathrm{i}\mathrm{x}\left(g\right)\right|$. A colouring $w$ is fixed by $g$ iff it is constant on each cycle of $g$’s action on $\Gamma$; the number of such colourings is $q^{\hspace{0.17em}c\left(g\right)}$, one free colour choice per cycle. Summing gives the stated formula, and rewriting the cycle-count distribution as a cycle-index polynomial yields the evaluated cycle index. $◼$

The building blocks of $c\left(g\right)$ are elementary and let the practitioner evaluate the sum without enumerating all bins. On a single axis of length $m$, a non-reflected axis contributes the identity permutation ($m$ fixed points) while a reflected axis contributes $\lfloor m/2\rfloor$ transpositions and $m\hspace{0.25em}\mathrm{m}\mathrm{o}\mathrm{d}\hspace{0.25em}2$ fixed points; the cell-cycle structure of a signed axis permutation is the product of these one-dimensional structures composed with the axis permutation $\pi$, exactly as in the wreath-product cycle index of $B_d$. We confirmed the one-dimensional reflection count for $m=2,\dots ,6$ and, more importantly, verified that the Burnside formula agrees with brute-force orbit counting in every small case: for instance $N\left(2,3,2\right)=102$, $N\left(2,2,3\right)=21$, $N\left(3,2,2\right)=22$, and $N\left(2,3,3\right)=2862$, each matching a direct orbit enumeration. Table 5 lists representative values, and Figure 4 contrasts the orbit count with the raw count $q^{\hspace{0.17em}{m}^d}$, showing the expected reduction by a factor approaching $\left|B_d\right|$ as the configurations become generic.

The cycle index via signed cycle types. The Burnside sum of Theorem 3 is governed by the cycle count $c\left(g\right)$, and $c\left(g\right)$ depends only on the signed cycle type of $g$: the multiset of pairs $\left(\mathcal{l},s\right)$, one per cycle of the axis permutation $\pi$, where $\mathcal{l}$ is the cycle length and $s\in \{+1,-1\}$ is the product of the reflection signs around that cycle. That the induced cycle structure of a signed permutation is determined by this data is classical: the cycle index of $B_d$ was first obtained by Harrison and High [33] and given in a simpler form by Chen [34], who computes the induced cycle structure on the vertex set of the hypercube $Q_d=\{0,1{\}}^d$ from the signed cycle decomposition. What we need here is the same computation for the grid ${\left[m\right]}^d$ with arbitrary resolution $m$; the case $m=2$ recovers the classical hypercube-vertex setting. We record the specialisation for completeness, since the $m$-dependence is what drives the enumeration below.

Lemma 1 (the block map and its cycle-length multiset). Let a cycle of $\pi$ have length $\mathcal{l}$ and net sign $s$ . On the corresponding factor ${\left[m\right]}^{\mathcal{l}}$ it induces the map

$$T\left(a_0,\dots ,a_{\mathcal{l}-1}\right)=\left(b_0,\dots ,b_{\mathcal{l}-1}\right), b_{\left(i+1\right)\hspace{0.25em}\mathrm{m}\mathrm{o}\mathrm{d}\hspace{0.25em}\mathcal{l}}=a_i, \mathrm{followed}\mathrm{by}{b}_0\mapsto \left(m-1\right)-b_0\hspace{0.25em}\mathrm{if}s=-1.$$

Then$T$is a permutation of${\left[m\right]}^{\mathcal{l}}$of order$\mathcal{l}$when$s=+1$and$2\mathcal{l}$when$s=-1$; write$L\left(\mathcal{l},s,m\right)$for its multiset of cycle lengths. Every cycle length divides the order of$T$, and the multiplicities are recovered from the fixed-point counts by Möbius inversion,

$$\#\left\{\mathrm{cycles}\mathrm{of}\mathrm{length}\mathrm{exactly}k\right\}={\displaystyle \frac{1}{k}}\sum _{e\mid k}\mu \left(k/e\right)\hspace{0.17em}\left|\mathrm{F}\mathrm{i}\mathrm{x}\left(T^e\right)\right|.$$

For $s=+1$ these counts are $\left|\mathrm{F}\mathrm{i}\mathrm{x}\left(T^k\right)\right|=m^{\mathrm{g}\mathrm{c}\mathrm{d}\left(k,\mathcal{l}\right)}$ , so the total number of cycles is $1/\mathcal{l}\sum _{k\mid \mathcal{l}}\phi \left(k\right)\hspace{0.17em}{m}^{\mathcal{l}/k}$ . For $s=-1$ the wrap-edge reflection makes the counts depend on the parity of $m$ , and $L\left(\mathcal{l},-1,m\right)$ is obtained from the same inversion. Explicitly, $T^k$ fixes a point exactly when the point is constant along the orbits of the $k$ -fold shift and, when $k/gcd(k,\mathcal{l})$ is odd, additionally reflection-symmetric in each such orbit, so $\left|Fix\left(T^k\right)\right|=m^{gcd(k,\mathcal{l})}$ when $k/gcd(k,\mathcal{l})$ is even, and ${\left(m\mathrm{m}\mathrm{o}\mathrm{d}2\right)}^{gcd(k,\mathcal{l})}$ otherwise; the cycle multiplicities of $L\left(\mathcal{l},-1,m\right)$ then follow by the same Möbius inversion. For $\mathcal{l}=1$ this gives the familiar pair $\left|Fix\left(T\right)\right|=m\mathrm{m}\mathrm{o}\mathrm{d}2$ and $\left|Fix\left(T^2\right)\right|=m$ . 

Proof of Lemma 1. $T$ is a composition of a slot rotation with (for $s=-1$) an involution of one coordinate, hence a permutation of the finite set ${\left[m\right]}^{\mathcal{l}}$. For $s=+1$, $T^{\mathcal{l}}$ is the identity and no smaller power is, so $\mathrm{o}\mathrm{r}\mathrm{d}\left(T\right)=\mathcal{l}$; a point is fixed by $T^k$ exactly when it is constant along the orbits of the shift by $k$, of which there are $\mathrm{g}\mathrm{c}\mathrm{d}\left(k,\mathcal{l}\right)$, giving $\left|\mathrm{F}\mathrm{i}\mathrm{x}\left(T^k\right)\right|=m^{\mathrm{g}\mathrm{c}\mathrm{d}\left(k,\mathcal{l}\right)}$, and the total cycle count is the classical necklace count $1/\mathcal{l}\sum _{k\mid \mathcal{l}}\phi \left(k\right)m^{\mathcal{l}/k}$. For $s=-1$, $T^{\mathcal{l}}$ applies the reflection $a\mapsto \left(m-1\right)-a$ to every coordinate—an involution with $m\hspace{0.25em}\mathrm{m}\mathrm{o}\mathrm{d}\hspace{0.25em}2$ fixed points—and $T^{2\mathcal{l}}=\mathrm{i}\mathrm{d}$, so $\mathrm{o}\mathrm{r}\mathrm{d}\left(T\right)=2\mathcal{l}$ and the fixed-point counts, hence the cycle lengths, depend on the parity of $m$. In both cases the cycle lengths divide $\mathrm{o}\mathrm{r}\mathrm{d}\left(T\right)$ and the stated Möbius inversion is the standard recovery of exact-period counts from fixed-point counts. $◼$

Theorem 4 (cycle index of$B_d$on the grid${\left[m\right]}^d$ ). For$g\in B_d$with signed cycle type$\{\left({\mathcal{l}}_1,s_1\right),\dots ,\left({\mathcal{l}}_r,s_r\right)\}$ , let$L_j=L\left({\mathcal{l}}_j,s_j,m\right)$be the cycle-length multiset of the block map of Lemma 1. Then the number of cycles of$g$on the full grid${\left[m\right]}^d$is

$$c\left(g\right)\hspace{0.33em}=\sum _{\left(a_1,\dots ,a_r\right)\in L_1\times \dots \times L_r}{\displaystyle \frac{a_1{a}_2\dots a_r}{\mathrm{l}\mathrm{c}\mathrm{m}\left(a_1,\dots ,a_r\right)}},$$

the product-action count taken over one cycle chosen from each block. Consequently$N\left(d,m,q\right)=1/2^{\hspace{0.17em}d}d!\sum _{g\in B_d}{q}^{c\left(g\right)}$is, for fixed$d$, a quasi-polynomial in$m$of degree$d$in the exponent with period dividing$2$. 

Proof. Distinct $\pi$-cycles act on disjoint coordinate blocks, so $g$ is the direct product of the block maps and the grid factorises as ${\left[m\right]}^d=\prod _j{\left[m\right]}^{{\mathcal{l}}_j}$. For a direct product of permutations, one cycle chosen from each factor with lengths $a_1,\dots ,a_r$ generates on their joint support exactly $\left(\prod _j{a}_j\right)/\mathrm{l}\mathrm{c}\mathrm{m}\left(a_1,\dots ,a_r\right)$ cycles—the orbit count of the diagonal generator of ${\mathbb{Z}}_{a_1}\times \dots \times {\mathbb{Z}}_{a_r}$—and summing over all tuples gives $c\left(g\right)$. The $s=+1$ block count is the necklace formula stated in Lemma 1; the $s=-1$ block introduces a reflection on the wrap edge whose effect depends on the parity of $m$, so the exponents are quasi-polynomials of period $2$, with leading term $m^d$ from the identity. $◼$

The formula is the evaluation at $a_k=q$ of the cycle-index polynomial $Z_{B_d}^{\Gamma }$ of the bin permutation representation. Read representation-theoretically, $q^{c\left(g\right)}$ counts the $g$-fixed colourings, so $N\left(d,m,q\right)$ is the dimension of the space of $B_d$-invariant functions on $A^{\Gamma }$—the invariant-theoretic quantity whose graded analogue is computed by a Molien series. We do not pursue the graded version, but it is the natural route to refinements such as counting occupancy classes by the number of occupied bins.

The formula specialises the classical hypercube computation of [33,34] to arbitrary $m$, and we verified it exhaustively: for all $g\in B_2$ and all $g\in B_3$ and every resolution $m\le 4$, the product-action formula reproduces the directly counted $c\left(g\right)$ in all $168$ cases.

Worked cycle index at $d=2$ . We derive $N\left(2,m,q\right)$ in closed form to illustrate the cycle-index computation. The group $B_2=D_4$ has eight elements acting on the $m^2$ bins: the identity $e$; the two axis reflections $r_x,r_y$; their product $r_x{r}_y$ (the ${180}^{\circ }$ rotation); the axis swap $s$ (matrix transpose); and the three remaining signed swaps $s{r}_x$, $s{r}_y$, $s{r}_x{r}_y$ (the $\pm {90}^{\circ }$ rotations and the anti-diagonal reflection). Counting cycles on the $m\times m$ board:

• $e$ fixes all $m^2$ bins, contributing $q^{m^2}$.
• Each axis reflection acts independently on each of the $m$ lines parallel to it, splitting the $m$ bins of every such line into $\lfloor m/2\rfloor$ transpositions together with $m\mathrm{m}\mathrm{o}\mathrm{d}2$ fixed bins, i.e. $\lceil m/2\rceil$ cycles per line, giving $c=m\lceil m/2\rceil$; the two such elements contribute $2q^{\hspace{0.17em}m\lceil m/2\rceil }$.
• The ${180}^{\circ }$ rotation fixes the centre bin when $m$ is odd and pairs the rest, giving $c=\lceil m^2/2\rceil$, contributing $q^{\lceil m^2/2\rceil }$.
• The transpose $s$ fixes the $m$ diagonal bins and pairs the $\left({\displaystyle \genfrac{}{}{0pt}{}{m}{2}}\right)$ off-diagonal transposed pairs, giving $c=m+\left({\displaystyle \genfrac{}{}{0pt}{}{m}{2}}\right)=\left({\displaystyle \genfrac{}{}{0pt}{}{m+1}{2}}\right)$; likewise the anti-diagonal reflection. These two contribute $2q^{\left({\displaystyle \genfrac{}{}{0pt}{}{m+1}{2}}\right)}$.
• The two ${90}^{\circ }$ rotations have $4$-cycles off the centre and one fixed centre when $m$ is odd, giving $c=\lceil m^2/4\rceil$; they contribute $2q^{\lceil m^2/4\rceil }$.

Hence

$$N\left(2,m,q\right)={\displaystyle \frac{1}{8}}\left(q^{m^2}+2\hspace{0.17em}{q}^{\hspace{0.17em}m\lceil m/2\rceil }+q^{\lceil m^2/2\rceil }+2\hspace{0.17em}{q}^{\left({\displaystyle \genfrac{}{}{0pt}{}{m+1}{2}}\right)}+2\hspace{0.17em}{q}^{\lceil m^2/4\rceil }\right).$$

For $m=3$, $q=2$ this evaluates to $1/8\left(2^9+2\cdot 2^6+2^5+2\cdot 2^6+2\cdot 2^3\right)=1/8\left(512+128+32+128+16\right)=102$, matching the exhaustive orbit count. For $m=3$, $q=3$ it gives $2862$, and for $m=4$, $q=2$ it gives $8548$; both agree with direct enumeration (Table 5). At q = 2 the formula counts binary square arrays up to the symmetries of the square; the resulting sequence is recorded as A054247 in the OEIS [37], and our expression reproduces every term listed there. We have not found a tabulation for general alphabet size, but we claim no novelty for the two-dimensional case: it is included as a worked illustration of the cycle-index computation.

A fully worked micro-example ( $d=2$ , $m=2$ ). The smallest non-trivial case can be checked entirely by hand and illustrates every ingredient. Here $\Gamma$ has $2^2=4$ bins, $B_2=D_4$ has $8$ elements, and with binary occupancy $q=2$ there are $2^4=16$ raw patterns. Burnside gives $N\left(2,2,2\right)=6$, and the six orbits are listed in Table 6 together with their sizes and stabiliser orders. The orbit sizes $1,4,4,2,4,1$ sum to $16$, and each satisfies the orbit–stabiliser relation $\left|\mathrm{o}\mathrm{r}\mathrm{b}\mathrm{i}\mathrm{t}\right|\cdot \left|\mathrm{s}\mathrm{t}\mathrm{a}\mathrm{b}\right|=8$: for example the two-bins-on-a-diagonal pattern has stabiliser of order $4$ (the diagonal is fixed setwise by the identity, the two reflections across the diagonals, and the ${180}^{\circ }$ rotation) and hence an orbit of size $2$.

The relation $\sum _{\mathrm{orbits}}|\mathrm{o}\mathrm{r}\mathrm{b}\mathrm{i}\mathrm{t}|=q^{m^d}$ visible in the table is a special case of the orbit-counting decomposition, and provides an independent arithmetic check on any Burnside computation.

Two further structural facts are worth recording. First, for fixed $m$ the count $N\left(d,m,q\right)$ is a polynomial in $q$ of degree $m^d$ whose coefficients are non-negative rationals summing appropriately; for instance $N\left(2,2,q\right)=1/8\left(q^4+2q^3+3q^2+2q\right)$ and $N\left(2,3,q\right)=1/8\left(q^9+4q^6+q^5+2q^3\right)$, which is the generating function of occupancy classes graded by alphabet size. Second, the orbit sizes are not uniform: at $d=2$, $m=2$, $q=2$ the six orbits split as two of size $1$, one of size $2$, and three of size $4$, whereas by $m=3$ the orbit-size distribution is $8$ orbits of size $1$, $8$ of size $2$, $50$ of size $4$ and $36$ of size $8$, so $66$ of the $102$ orbits still have non-trivial stabiliser—a proportion that, by Proposition 2, tends to zero as $q$ or $m$ grows. The small grids are thus the exception in which symmetry is abundant; the generic large-grid instance sits in a full-size orbit.

Figure 5 shows these distributions directly.

Proposition 1 (orbit sizes divide the group order). For every occupancy function $w$ , the orbit $B_d\cdot w$ has size $\left|B_d\right|/\left|\mathrm{S}\mathrm{t}\mathrm{a}\mathrm{b}\left(w\right)\right|$ , a divisor of $2^{\hspace{0.17em}d}d!$ ; in particular every orbit size divides $2^{\hspace{0.17em}d}d!$ , and the generic orbit (trivial stabiliser) has the full size $2^{\hspace{0.17em}d}d!$ . 

Proof. This is the orbit–stabiliser theorem applied to $B_d$ acting on $A^{\Gamma }$: the map $g\mathrm{S}\mathrm{t}\mathrm{a}\mathrm{b}\left(w\right)\mapsto g\cdot w$ is a bijection from left cosets of $\mathrm{S}\mathrm{t}\mathrm{a}\mathrm{b}\left(w\right)$ onto the orbit, so $|B_d\cdot w|=\left[B_d:\mathrm{S}\mathrm{t}\mathrm{a}\mathrm{b}\left(w\right)\right]=\left|B_d\right|/\left|\mathrm{S}\mathrm{t}\mathrm{a}\mathrm{b}\left(w\right)\right|$, which divides $\left|B_d\right|=2^{\hspace{0.17em}d}d!$ by Lagrange’s theorem. When $\mathrm{S}\mathrm{t}\mathrm{a}\mathrm{b}\left(w\right)$ is trivial the orbit attains the full size. $◼$

Corollary 2 (closed form at$d=3$ , the octahedral case).Applying Theorem 4 to$B_3=O_h$ —whose$48$elements fall into ten conjugacy classes (identity, face and edge reflections,$\pm {90}^{\circ }$and${180}^{\circ }$face rotations,${120}^{\circ }$vertex rotations, the central inversion, and their compositions)—yields an explicit closed form for every resolution, split by the parity of$m$ . For$m$even,

$$N\left(3,m,q\right)={\displaystyle \frac{1}{48}}\left(q^{m^3}+13\hspace{0.17em}{q}^{m^3/2}+6\hspace{0.17em}{q}^{m^3/2+m^2/2}+12\hspace{0.17em}{q}^{m^3/4}+8\hspace{0.17em}{q}^{m^3/3+2m/3}+8\hspace{0.17em}{q}^{m^3/6+m/3}\right),$$

and for $m$ odd,

$$N\left(3,m,q\right)={\displaystyle \frac{1}{48}}\left(q^{m^3}+9\hspace{0.17em}{q}^{m^3/2+m^2/2}+9\hspace{0.17em}{q}^{m^3/2+m/2}+q^{m^3/2+1/2}+6\hspace{0.17em}{q}^{m^3/4+3m/4}+6\hspace{0.17em}{q}^{m^3/4+m/4+1/2}+8\hspace{0.17em}{q}^{m^3/3+2m/3}+8\hspace{0.17em}{q}^{m^3/6+m/3+1/2}\right).$$

Setting $m=3$ gives $N\left(3,3,2\right)=2,852,288$ and $N\left(3,3,3\right)=158,942,078,604$; setting $m=2$ gives $N\left(3,2,2\right)=22$. All three, and the values up to $m=5$, $q=3$ (a $58$-digit integer), were confirmed against independent orbit enumeration. Two remarks place this formula. First, at $m=2$ the grid is the vertex set of the cube and the count reduces to the classical enumeration of $q$-colourings of the cube vertices under $O_h$, giving $N\left(3,2,2\right)=22$; that case is standard. Second, for $m\ge 3$—the regime that matters for clustering grids—we are not aware of a published closed form, and the parity split above appears not to have been recorded, though the ingredients are entirely classical. Unlike its two-dimensional counterpart, the sequence of values it produces returns no match in the OEIS [37] (searched 28 July 2026), whereas the two-dimensional sequence is recorded there as A054247. We therefore present it as an explicit specialisation rather than a new enumerative technique. In general, as the resolution $m$ or the alphabet $q$ grows, the identity term $q^{m^d}$ dominates the Burnside sum and

$$N\left(d,m,q\right)\hspace{0.33em}=\hspace{0.33em}{\displaystyle \frac{q^{m^d}}{2^{\hspace{0.17em}d}d!}}\hspace{0.17em}\left(1+o\left(1\right)\right),$$

so the number of distinct configurations approaches the raw count divided by the full group order $\left|B_d\right|=2^{\hspace{0.17em}d}d!$. We state this precisely.

Proposition 2 (asymptotic orbit count). For fixed $d$ , as $q\to \mathrm{\infty }$ (or as $m\to \mathrm{\infty }$ with $q\ge 2$ ), $N\left(d,m,q\right)={\displaystyle \frac{q^{m^d}}{2^{\hspace{0.17em}d}d!}}\left(1+o\left(1\right)\right)$ . 

Proof. The identity element contributes $q^{m^d}$ to the Burnside sum, while every non-identity $g$ has $c\left(g\right)\le m^d-1$ (it moves at least one bin, merging at least two of them into a common cycle). More sharply, every non-identity element satisfies $c\left(g\right)\le c_2\left(d,m\right)$ with $c_2\left(d,m\right)$ the largest non-identity cycle count identified below, so the total non-identity contribution is at most $\left(2^{\hspace{0.17em}d}d!-1\right)\hspace{0.17em}{q}^{c_2\left(d,m\right)}$, which is an $O\left(1/q\right)$ fraction of the identity term whenever $m^d-c_2\left(d,m\right)\ge 1$. Both limiting regimes now follow. For fixed $m$, the ratio is $O\left(1/q\right)$ and vanishes as $q\to \mathrm{\infty }$. For fixed $q\ge 2$, the gap $m^d-c_2\left(d,m\right)$ is itself unbounded in $m$: the largest non-identity cycle count is attained by a single axis transposition, which fixes the $m^{d-1}$ bins whose two transposed coordinates agree and pairs up the rest, giving $c_2\left(d,m\right)=\left(m^d+m^{d-1}\right)/2$ and hence a gap of $\left(m^d-m^{d-1}\right)/2$, which tends to infinity with $m$ for every fixed $d\ge 2$. Consequently $q^{-(m^d-c_2)}$ tends to zero as well. Dividing by $\left|B_d\right|$ gives the claim.   ■

Remark 1 (the second-order term). The correction can be made completely explicit. Write $c_2\left(d,m\right)$ for the largest cycle count attained by a non-identity element and $\mu \left(d,m\right)$ for the number of elements attaining it. Then the two-term expansion reads

$$N\left(d,m,q\right)\hspace{0.33em}=\hspace{0.33em}{\displaystyle \frac{1}{2^{\hspace{0.17em}d}d!}}\left(q^{m^d}+\mu \hspace{0.17em}{q}^{\hspace{0.17em}{c}_2\left(d,m\right)}+O\left(q^{\hspace{0.17em}{c}_3\left(d,m\right)}\right)\right),$$

and the relative error is therefore

$$\left|\hspace{0.17em}{\displaystyle \frac{N\left(d,m,q\right)}{q^{m^d}/\left|B_d\right|}}-1\hspace{0.17em}\right|\hspace{0.33em}=\hspace{0.33em}\mu \left(d,m\right)\hspace{0.17em}{q}^{-\left(m^d-c_2\left(d,m\right)\right)}\left(1+O\left(1/q\right)\right).$$

At $d=2$ we have $c_2\left(2,m\right)=\left({\displaystyle \genfrac{}{}{0pt}{}{m+1}{2}}\right)$, so the decay exponent is $m^2-\left({\displaystyle \genfrac{}{}{0pt}{}{m+1}{2}}\right)=\left({\displaystyle \genfrac{}{}{0pt}{}{m}{2}}\right)$; the value is attained by the two diagonal reflections for every $m$, joined by the two axis reflections when $m$ is odd, giving $\mu =2$ for $m$ even and $\mu =4$ for $m$ odd. At $d=3$ we have $c_2\left(3,m\right)=1/2\left(m^3+m^2\right)$ and decay exponent $1/2\left(m^3-m^2\right)$; here the value is attained by the six diagonal mirror planes—the three pure axis transpositions together with the three that also flip the two swapped signs—for every $m$, joined by the three face reflections when $m$ is odd, giving $\mu =6$ for $m$ even and $\mu =9$ for $m$ odd. Numerically, at $d=2$, $m=5$, $q=8$ the neglected terms contribute a relative $3.7\times {10}^{-9}$ against the predicted $4\cdot 8^{-10}=3.7\times {10}^{-9}$, and at $d=3$, $m=5$, $q=8$ a predicted $9\cdot 8^{-50}=6.3\times {10}^{-45}$. $◼$

This asymptotic is the precise statement that “almost all” occupancy functions have trivial stabiliser and therefore full-size orbits—the same fact that, read from the algorithmic side, tells us there is no per-instance bin collapse to exploit on generic data (Section 7). Figure 4 shows the finite-$m$ approach to this ratio.

The case $d=2$, $m=3$, $q=2$ is a useful anchor: there are $102$ distinct empty/non-empty bin patterns on a $3\times 3$ grid up to $B_2=D_4$. This is the clustering-grid counterpart of the DR Code’s orbit counts on the same $3\times 3$ board, with the recoverability constraint removed and the alphabet reduced to “empty vs. occupied”.

7. Canonical Form and Complexity

Enumeration counts orbits; applications need a way to name them. We give a canonical form—a function constant on orbits and injective across them—together with its cost.

Definition 3 (canonical form). Order $A^{\Gamma }$ lexicographically by reading bin values in a fixed bin order. For an occupancy function $w$ define

$$\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w\right)\hspace{0.33em}=\hspace{0.33em}\underset{g\in B_d}{\mathrm{m}\mathrm{i}\mathrm{n}}\hspace{0.33em}\left(g\cdot w\right),$$

the lexicographically smallest function in the orbit of $w$, and let $\mathrm{r}\mathrm{e}\mathrm{p}\left(w\right)\in B_d$ be a group element attaining the minimum.

Algorithm 1 makes the computation in Definition 3 explicit. It enumerates the $\left|B_d\right|$ signed permutations, applies each to the occupancy function, and keeps the lexicographically smallest image together with the group element that produced it.

Algorithm 1: Canonical form of an occupancy function under $B_d$

Input: occupancy function $w:\Gamma \to A$on the grid $\Gamma ={\left[m\right]}^d$

Output: canonical form $w^{\star }=\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w\right)$and a witness $g^{\star }\in B_d$with $g^{\star }\cdot w=w^{\star }$

1. $w^{\star }\leftarrow \perp$; $g^{\star }\leftarrow \perp$

2. for each signed permutation $g=\left(\pi ,\epsilon \right)\in B_d$ do

3.    compute the push-forward $u\leftarrow g\cdot w$(relabel all $m^d$bins)

4.    if $w^{\star }=\perp$ or $u{<}_{\mathrm{l}\mathrm{e}\mathrm{x}}{w}^{\star }$ then $w^{\star }\leftarrow u$; $g^{\star }\leftarrow g$

5. return $\left(w^{\star },g^{\star }\right)$

The loop runs exactly $\left|B_d\right|=2^{\hspace{0.17em}d}d!$ times; each iteration relabels $m^d$ bins and performs one lexicographic comparison of length $m^d$, so the total cost is $O\left(2^{\hspace{0.17em}d}d!\hspace{0.17em}{m}^d\right)$ time and $O\left(m^d\right)$ working space, as stated in Theorem 5. Because the group is enumerated in full, the algorithm is deterministic and its worst case coincides with its every case; there is no data-dependent search. Two complementary remarks bound this cost from below and from above in the typical case. For the lower bound, any correct canonical form must inspect essentially the whole grid: at $d=2$, $m=3$, $q=2$, flipping any single one of the nine bins changes the canonical form; the same holds on any grid, since altering one bin changes the multiset of occupancy values, which is an orbit invariant. Consequently no algorithm can decide the canonical representative after reading $o\left(m^d\right)$ bins, and the $m^d$ factor in Theorem 5 is unavoidable—only the group factor $2^{\hspace{0.17em}d}d!$ is specific to our method. For the typical case, the lexicographic comparison in step 4 aborts at the first differing bin, which on random occupancy data happens after $O\left(1\right)$ bins; measured over random instances the total work falls to a little over half of the naive bound that charges a full-length comparison to every group element (0.51–0.64 across the parameter range we tested), with the $\left|B_d\right|$ push-forwards, not the comparisons, dominating.

Theorem 5 (canonical form: correctness and complexity).The map$\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}$satisfies$\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(g\cdot w\right)=\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w\right)$for all$g\in B_d$, and$\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w\right)=\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w\mathrm{\text{'}}\right)$if and only if$w\mathrm{\text{'}}$lies in the orbit of$w$. It can be computed in$O\left(2^{\hspace{0.17em}d}d!\hspace{0.17em}\cdot m^{\hspace{0.17em}d}\right)$time and$O\left(m^{\hspace{0.17em}d}\right)$working space by evaluating the$\left|B_d\right|$push-forwards and retaining the lexicographic minimum. The number of distinct canonical forms over all occupancy functions equals$N\left(d,m,q\right)$of Theorem 3. 

Proof. For invariance, the orbit of $g\cdot w$ equals the orbit of $w$, so the minimum over the orbit is the same function; hence $\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(g\cdot w\right)=\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w\right)$. For separation, if $\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w\right)=\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w\mathrm{\text{'}}\right)$ then $w$ and $w\mathrm{\text{'}}$ share an orbit element and therefore share an orbit; conversely if they share an orbit they share its unique lexicographic minimum. The algorithm iterates over the $\left|B_d\right|=2^{\hspace{0.17em}d}d!$ group elements; each push-forward relabels $m^{\hspace{0.17em}d}$ bins and the lexicographic comparison is linear in $m^{\hspace{0.17em}d}$, giving $O\left(2^{\hspace{0.17em}d}d!\hspace{0.17em}{m}^{\hspace{0.17em}d}\right)$ time, while only the current best function and the candidate are stored, giving $O\left(m^{\hspace{0.17em}d}\right)$ space. Since $\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}$ is constant on orbits and distinct across them, the number of images equals the number of orbits, which is $N\left(d,m,q\right)$ by Theorem 3. $◼$

Two consequences are worth stating separately, since they are what a downstream system actually uses.

Corollary 3 (uniqueness of the normal form).Each$B_d$-orbit contains exactly one occupancy function that is its own canonical form, namely$w^{\star }=\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w\right)$for any$w$in the orbit; and$\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w^{\star }\right)=w^{\star }$. Hence$\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}$is a well-defined normal form: it selects one and only one representative per orbit, and applying it twice changes nothing. 

Proof. Existence: for $w$ in the orbit, $w^{\star }=\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w\right)$ lies in the orbit (it is $g\cdot w$ for some $g\in B_d$) and is the lexicographic minimum of that orbit. Idempotence: since $w^{\star }$ lies in the same orbit, $\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w^{\star }\right)$ is the minimum of the same finite totally ordered set, so $\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w^{\star }\right)=w^{\star }$. Uniqueness: if $u$ in the orbit satisfies $\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(u\right)=u$ then $u$ is the minimum of the orbit, and a totally ordered finite set has exactly one minimum, so $u=w^{\star }$. $◼$

Corollary 4 (orbit membership test).For occupancy functions$w,w\mathrm{\text{'}}$on the same grid,$w\mathrm{\text{'}}\in B_d\cdot w$if and only if$\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w\right)=\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w\mathrm{\text{'}}\right)$ . The test therefore costs$O\left(2^{\hspace{0.17em}d}d!\hspace{0.17em}{m}^d\right)$time in general, and only$O\left(m^d\right)$ —a single key comparison—when canonical forms have already been computed and stored, as they are in the indexing and caching applications ofSection 8.

Proof. If $w\mathrm{\text{'}}=g\cdot w$ then the two orbits coincide and hence so do their lexicographic minima. Conversely, equal canonical forms mean $w$ and $w\mathrm{\text{'}}$ share an orbit member, and orbits are either equal or disjoint, so the orbits coincide and $w\mathrm{\text{'}}\in B_d\cdot w$. The cost is that of two applications of Algorithm 1 plus one comparison of length $m^d$, which is $O\left(2^{\hspace{0.17em}d}d!\hspace{0.17em}{m}^d\right)$; when the canonical forms are cached, only the comparison remains. $◼$

We verified Corollary 4 directly: over $600$ random pairs at $d\in \{2,3\}$, the canonical-form test agreed with an exhaustive search over all $g\in B_d$ in every case.

Two features distinguish this canonical form from general graph or molecule canonicalisation: the group is fixed and explicit, so the loop length $2^{\hspace{0.17em}d}d!$ depends only on the ambient dimension and no search over an unknown automorphism group is needed; and the cost is modest in the low dimensions where grid clustering is used—$\left|B_2\right|=8$, $\left|B_3\right|=48$, $\left|B_4\right|=384$. Table 7 lists the group size and resulting cost in the low dimensions relevant to grid clustering. We verified correctness computationally—exhaustively over all binary occupancy functions of the smallest grid, and on random instances for $d\in \{2,3\}$: applying any $g$ before canonicalising never changed the result, and functions sharing a canonical form always shared an orbit.

A worked $3\times 3$ example. Consider $d=2$, $m=3$, binary occupancy $q=2$ (empty/occupied). There are $2^9=512$ raw occupancy patterns and, by the formula above, $102$ orbits under $B_2=D_4$. Take the pattern with occupied bins $\{\left(0,0\right),\left(0,1\right),\left(1,0\right)\}$ (an L-shaped corner cluster). Its $D_4$-orbit consists of this corner L pointed into each of the four corners of the board, four distinct patterns in all: the pattern is fixed by the reflection across the main diagonal through its occupied corner, which swaps $\left(0,1\right)\leftrightarrow \left(1,0\right)$ and fixes $\left(0,0\right)$, so its stabiliser has order $2$ and its orbit has size $8/2=4$ by the orbit–stabiliser theorem. The canonical representative is the lexicographically smallest of the four; because empty bins carry the smallest alphabet value, the minimum is the image that pushes the occupied bins as late as possible in the bin order, namely the L occupying the opposite corner, and every one of the four input patterns canonicalises to it. A downstream index that keys clustering results by this canonical corner stores and computes the result once rather than four times. A generic (asymmetric) occupancy pattern would instead have trivial stabiliser and a full orbit of size $8$; the corner L is one of the $66$ of the $102$ orbits with non-trivial stabiliser recorded in Section 6, which is precisely why its orbit is smaller. This is the mechanism of Theorem 5 in miniature, and it is exactly what a $B_d$-closed augmentation of a benchmark suite triggers at scale.

The three-dimensional case is equally concrete. On the $2\times 2\times 2$ cube ($d=3$, $m=2$) under $B_3=O_h$, a single occupied corner has orbit size $8$ (the eight cube corners) and stabiliser of order $6$—the subgroup $S_3$ permuting the three axes that meet at that corner. Two occupied corners sharing an edge form an orbit of size $12$ with stabiliser of order $4$, and a fully occupied face (four corners) has orbit size $6$ with stabiliser of order $8$. These are exactly the vertex, edge, and face orbits of the cube, recovered here as occupancy classes; canonicalising any member of each orbit returns the same lexicographically minimal representative.

Scope of the benefit. The canonical form does not accelerate a single clustering run: on a generic occupancy function the stabiliser is trivial, the orbit has full size $2^{\hspace{0.17em}d}d!$, and there is nothing to collapse within one instance. We measured the fraction of random occupancy functions with a non-trivial stabiliser and found it falls rapidly toward zero as the grid and alphabet grow: at $d=3$, $m=3$, $q=4$ it was already $0\mathrm{\%}$ over the sampled instances (Figure 6a). The orbit structure pays off, instead, across collections of instances. In a benchmark library closed under $B_d$—for example a synthetic suite augmented by axis permutations and reflections, a common way to test rotation/reflection invariance—canonicalising removes the redundant members: a library of $400$ instances built as full orbits of $50$ seeds collapses to exactly $50$ canonical instances, an $8.00\times$ reduction equal to the average orbit size, drawn as the dashed reference line in Figure 6b. That factor is a construction, not an observation: the library was built from complete orbits, so it verifies that the bound of Definition 4 is attained and nothing more. The complementary question—whether symmetry-equivalent instances arise by chance in a workload that was not built that way—we settled by measurement, and the answer is essentially no. Over 20,000 independently drawn occupancy functions we compared the number of distinct functions with the number of distinct canonical forms, so as to isolate what canonicalisation adds beyond hashing the occupancy function itself. The additional factor was 5.02 at d = 2, m = 3, q = 2 and 2.25 at m = 4, but it fell to 1.002 at m = 5, 1.000 at m = 6, 1.004 at d = 3, m = 3, q = 2, and 1.000 at both d = 2, m = 4, q = 4 and d = 3, m = 3, q = 4 (Figure 6b). The early values are an artefact of a workload larger than the instance space: at d = 2, m = 3, q = 2 the 20,000 draws saturate all 512 patterns, which collapse to the 102 orbits of Table 5. Once the space is large relative to the workload—which is the case at every resolution of practical interest—symmetry-based deduplication adds nothing to ordinary hashing. It pays only on workloads that are deliberately symmetry-closed, or that re-present the same measurement under a different axis convention. The same mechanism lets a pipeline that meets the same measurement again under a different axis convention reuse a cached clustering keyed by canonical form rather than recompute it. These are deduplication, indexing, and caching gains, and we claim only them.

8. Discussion

What the symmetry framework provides. The theorems combine into one practical statement: because clustering is equivariant (Theorem 2), instances equal up to axis relabelling are interchangeable, so replacing each by its canonical form (Theorem 5) is lossless, that representative is unique per orbit (Corollary 3), and the number of canonical instances is known in advance (Theorem 3). A pipeline can therefore store, index, and cache instances by canonical form with no effect on the result beyond a fixed relabelling of axes. Table 8 summarises where the framework helps and where it does not.

Figure 7 places the construction in a concrete workflow.

Definition 4 (symmetry-closed library). A finite collection $\mathcal{L}$ of occupancy functions is $B_d$ -closed if $g\cdot w\in \mathcal{L}$ for all $w\in \mathcal{L}$, $g\in B_d$; equivalently $\mathcal{L}$ is a union of orbits. Its canonical reduction $\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(\mathcal{L}\right)=\{\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(w\right):w\in \mathcal{L}\}$ holds exactly one member per orbit, so the deduplication factor $\left|\mathcal{L}\right|/\left|\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(\mathcal{L}\right)\right|$ equals the average orbit size of $\mathcal{L}$ and is at most $2^{\hspace{0.17em}d}d!$.

Algorithm 2 states the deduplication and caching procedure that Definition 4 supports.

Algorithm 2: Canonical deduplication and result caching

Input: a stream or library $\mathcal{L}$ of occupancy functions; a clustering routine $\mathcal{C}$

Output: the clustering of every $w\in \mathcal{L}$, with $\mathcal{C}$ evaluated once per orbit

1. initialise an empty table $T$ keyed by canonical forms

2. for each $w\in \mathcal{L}$ do

3.   $\left(w^{\star },g^{\star }\right)\leftarrow$ Algorithm 1 applied to $w$

4.    if $w^{\star }\notin T$ then $T\left[w^{\star }\right]\leftarrow \mathcal{C}\left(w^{\star }\right)$

5.    report $\mathcal{C}\left(w\right)={\left(g^{\star }\right)}^{-1}\cdot T\left[w^{\star }\right]$

6. return $T$

Step 5 is exactly Theorem 2, and Corollary 1 guarantees that any size-based summary of the cached clustering is already the correct value for $w$: the stored clustering of the canonical form is transported back to $w$ by the inverse of the witness element, with no recomputation. The number of calls to $\mathcal{C}$ equals the number of distinct orbits met, so on a $B_d$-closed library the saving is the average orbit size, and on a library with no symmetric repeats it is nil—the guarantee is one-sided. The table $T$ stores one key of $m^d$ symbols and one clustering per orbit, so its memory is $\left|\mathrm{c}\mathrm{a}\mathrm{n}\mathrm{o}\mathrm{n}\left(\mathcal{L}\right)\right|\cdot O\left(m^d\right)$, a factor of the average orbit size smaller than storing results per instance.

Benchmark deduplication. Synthetic suites are frequently built to test invariance by including, for each base dataset, its images under axis permutations and reflections; lifted to the bin grid, such a suite is closed under (a subgroup of) $B_d$. Indexing it by canonical form collapses each orbit to one representative, reducing storage and clustering runs by the average orbit size, up to $\left|B_d\right|$. Our controlled run realised this exactly: $400$ instances built as full $B_2$-orbits of $50$ seeds reduced to $50$ canonical instances, an $8.00\times$ saving (the dashed line of Figure 6b). That suite was constructed to be orbit-closed, so the figure confirms that the bound is attained rather than measuring how often symmetry occurs unbidden; on independently drawn instances it does not (Section 7). Re-analysis under a different axis convention. A measurement re-analysed after its channels are exported in a different order, or with an axis orientation reversed by a different instrument or compensation convention, presents a grid that coincides with one already seen up to axis relabelling. Keying a result cache by canonical form lets the pipeline reuse a previously computed clustering whenever a new grid is $B_d$-equivalent to one already seen, converting repeated work into a single canonical computation plus an $O\left(2^{\hspace{0.17em}d}d!\hspace{0.17em}{m}^d\right)$ lookup key. Deciding whether a newly arriving grid has been seen before is exactly the orbit membership test of Corollary 4, which is a single key lookup once canonical forms are stored. Instance registries. When clustering results are archived for audit or comparison—common in regulated settings such as clinical flow cytometry—the canonical form provides a stable, presentation-independent key under which equivalent instances are recognised as identical, preventing the same analysis from being logged as several distinct records merely because its axes were listed differently.

$\mathrm{up}\mathrm{to}\left|B_d\right|\times$In each scenario the canonical form is computed once per instance at cost $O\left(2^{\hspace{0.17em}d}d!\hspace{0.17em}{m}^d\right)$, a small constant multiple of the cost of reading the grid in the dimensions of practical interest. The saving is realised whenever instances recur up to symmetry and is zero when they do not, so the method never harms and helps exactly to the extent that symmetry is present in the workload. In practice the decision is easy to make in advance: apply canonicalisation whenever the workload is generated by augmentation, or by repeated presentation of the same measurement under differing axis conventions (where symmetric repeats are guaranteed), and skip it for one-off analyses of a single dataset (where Proposition 2 says there will be nothing to collapse). A break-even rule makes this quantitative: canonicalising one instance costs $O\left(2^{\hspace{0.17em}d}d!\hspace{0.17em}{m}^d\right)$, which is a small constant number of grid passes, so a cached pipeline becomes cheaper as soon as the number of symmetry-equivalent revisits exceeds that constant. In a $B_d$-closed library the deduplication factor equals the average orbit size, bounded by $2^{\hspace{0.17em}d}d!$; the controlled run reported in Section 7 attains the full factor $8.00$ at $d=2$.

Symmetry, partial symmetry, and symmetry breaking. It is worth stating plainly where the symmetry of this paper lives and how data destroys it. The grid $\Gamma ={\left[m\right]}^d$ is fully $B_d$-symmetric as a combinatorial object: every signed axis permutation is an automorphism of its adjacency graph. An occupancy function $w$ placed on that grid generally breaks this symmetry, and its stabiliser $\mathrm{S}\mathrm{t}\mathrm{a}\mathrm{b}\left(w\right)\le B_d$ measures exactly how much survives. The possible residual symmetries are therefore the subgroups of $B_d$, and the stabiliser is the precise invariant recording which one a given dataset realises. At $d=2$ the lattice is small enough to write down: $B_2=D_4$ has ten subgroups—the trivial group, five of order $2$ (four reflections and the ${180}^{\circ }$ rotation), three of order $4$ (the cyclic rotation group and two Klein four-groups), and $D_4$ itself—so any two-dimensional occupancy pattern falls into one of at most ten symmetry types, and the orbit-size distribution reported above is a coarsening of that classification. Not every subgroup is realisable on every grid: on the $3\times 3$ binary grid only nine of the ten occur, the cyclic rotation group $C_4$ being unattainable as an exact stabiliser because any pattern invariant under the ${90}^{\circ }$ rotation is constant on the corner, edge and centre orbits and hence also reflection-invariant; all ten occur once $m\ge 4$.

Figure 8 displays this lattice.

This is symmetry breaking in the ordinary group-theoretic sense, and it is the reason our claims are stated the way they are. Fully symmetric instances (stabiliser $B_d$, such as the empty or uniformly full grid) are combinatorially exceptional; partially symmetric instances (stabiliser a proper non-trivial subgroup) are common only on small grids; and generic instances break the symmetry completely, giving $\mathrm{S}\mathrm{t}\mathrm{a}\mathrm{b}\left(w\right)=\left\{e\right\}$ and a full-size orbit. Proposition 2 quantifies the last statement asymptotically. The practical consequence is the asymmetry noted in the abstract: the symmetry of the grid is what makes enumeration and canonical naming possible, while the asymmetry of the data is what removes any per-instance shortcut. Both halves are needed to state the contribution honestly, and both are properties of the same group action.

Limitations. Three boundaries delimit the contribution. First and most important, the benefit is across instances, not within one: a single generic dataset has trivial stabiliser, so canonicalisation neither reduces its bins nor speeds its clustering, and we make no such claim. Second, the analysis concerns the grid of bins, not the raw points; methods whose clusters depend on sub-bin geometry beyond bin occupancy and bin adjacency fall outside the model, though the large grid-based family we target does not. Third, the group is the pure geometric grid symmetry $B_d$; if an application also identifies feature axes by an external relation (a domain symmetry beyond axis relabelling), the relevant group is a supergroup of $B_d$ and the enumeration and canonical form must be recomputed for it, exactly as the DR Code augments $D_4$ by its combinatorial $C_3\times C_3$ factor.

Verification methodology. Every quantitative assertion in this paper was checked computationally, exhaustively wherever the relevant finite structure is small enough to enumerate and by large randomised trials otherwise. Theorem 1 was verified for $d\in \{2,3,4\}$ and $m\le 4$ (bijectivity and adjacency preservation on every pair; the group-composition identity was verified separately over all element pairs); Theorem 2 on $6,000$ random (occupancy, group-element, threshold) triples—$d=2$ with $m\in \{3,4,5\}$ and $d=3$ with $m\in \{3,4\}$; Theorem 3 by evaluating the Burnside sum against an independent brute-force orbit enumeration for six parameter triples small enough to enumerate directly, namely $\left(2,2,2\right)$, $\left(2,2,3\right)$, $\left(2,3,2\right)$, $\left(2,3,3\right)$, $\left(2,4,2\right)$ and $\left(3,2,2\right)$, with exact agreement (e.g. $N\left(2,3,2\right)=102$, $N\left(3,2,2\right)=22$); and Theorem 5 by confirming that pre-applying a group element never changes the canonical form and that two functions share a canonical form exactly when they share an orbit. The trivial-stabiliser limitation was quantified by sampling (Figure 6a), and the deduplication that canonicalisation adds over plain hashing was measured on 20,000 independently drawn occupancy functions across eight parameter settings (Figure 6b). These are computational certifications of finite mathematical statements rather than experiments on data. All checks ran in Python 3.12 using only the standard library, complete in under a minute on a commodity laptop, and are reproduced by the open-source script (see Data Availability); the randomised checks are large enough to be seed-independent and the enumeration checks are exhaustive.

Relation to the DR Code programme. This paper completes a natural arc. The 2D DR Code identified $D_4$ as the geometric symmetry of a constrained $3\times 3$ grid; the 3D DR Code raises it to $O_h$ on a $3\times 3\times 3$ cube [14]; the present work identifies both as $B_2$ and $B_3$ and treats the general $B_d$ acting on an unconstrained clustering grid. The methodological through-line—orbit decomposition, Burnside enumeration, canonical-form deduplication—is shared, but the object is different: an unbounded space of occupancy functions rather than a finite catalogue of recoverable layouts, and a purely geometric group rather than one fused with a code-specific constraint. The DR Code constraint itself places its layouts among Latin squares [38] and, more generally, frequency (F-)squares [39]; the clustering grid carries no such row/column constraint, which is precisely why the relevant group here is the pure geometric symmetry $B_d$.

9. Conclusions

We have given a symmetry analysis of grid-based density clustering built on the hyperoctahedral group $B_d=C_2\wr S_d$ of order $2^{\hspace{0.17em}d}d!$. We proved that $B_d$ acts on bin-occupancy functions by adjacency-preserving bijections (Theorem 1), that grid density clustering is $B_d$-equivariant (Theorem 2), that bin-occupancy configurations up to $B_d$ are counted by a Burnside average, for which we specialised the classical cycle index of $B_d$ to the grid and gave explicit closed forms in two and three dimensions—the three-dimensional formula being a parity-split expression, valid at every resolution, that we have not found recorded elsewhere (Theorems 3 and 4, Corollary 2), and that a fixed-group canonical form computes orbit representatives in $O\left(2^{\hspace{0.17em}d}d!\hspace{0.17em}{m}^{\hspace{0.17em}d}\right)$ time with correctness guaranteed by construction (Theorem 5). We were explicit that the payoff is canonical indexing, deduplication of symmetry-closed libraries, and cross-run caching, not a per-instance speedup, and we quantified each. The framework specialises to the $D_4$ and $O_h$ symmetries of the 2D and 3D DR Code, of which it is the unconstrained generalisation. Four extensions follow naturally. When an application identifies feature axes through an external relation, the relevant symmetry is a supergroup of $B_d$ and Section 6 and Section 7 carry over once that group is substituted. The occupancy alphabet can be made weighted or continuous; equivariance is unaffected, but the enumeration must then be replaced by an integral or by orbit counting on a quotient—an open combinatorial question. For very high-dimensional grids an approximate canonical form based on hashed bin representations would extend the method to regimes where $2^{\hspace{0.17em}d}d!$ is itself large. Finally, a principled canonical choice of the density threshold $\tau$ would make a clustering pipeline presentation-independent end to end. Beyond the specific setting, the broader point is a reproducibility one. Automated gating and other grid-clustering pipelines are increasingly used in regulated and clinical contexts, where an analysis is expected to be reproducible from its recorded inputs; an analysis whose output depends on the incidental order in which channels were written is not. Canonical axis handling removes that dependence at negligible cost and makes symmetry-equivalent analyses recognisable as the same analysis. All results reported here were confirmed computationally—by exhaustive enumeration on the relevant finite structures where that is feasible, and by large randomised trials elsewhere.