study the degree of aggregation as a function of $\beta$ and $\xi$. Whenever the system forms a steady-state cluster distribution, we explore the clusters morphology by computing their mean size, the cluster-size distribution, their radius of gyration, and both the polar and nematic order parameter. When clusters are present, we identify them following a distance criterion: at each time step, two particles belong to the same cluster of size $s$ whenever their distance is smaller than $r_{cl} = 1.75 \sigma$ (first minimum of the $g(r)$, see appendix \[Appendix:Clusters\]).

To calculate the cluster-size distribution $f(s)$ we apply a criterion based on [@chantalbins2012]: (i) We arbitrarily subdivide the range of $s$-values into intervals $\Delta s_i = (s_{i,max} - s_{i,min})$, where $n_i^t$ is the total number of clusters within each interval $\Delta s$; (ii) we assign the value $n_i=n^t_i/\Delta s_i$ to every $s$ within $\Delta s_i$, and compute the fraction of clusters of size $s$ as $f(s)=n_i/N_c$, where $N_c=\sum_i n_i \Delta s_i$ is the total number of clusters.

To characterise the clusters’ morphology, we compute the radius of gyration $$\label{eq:Rg}
R_g(s) = \sqrt{\sum_{i,j=1}^s \frac{\left(\vec{r}_i-\vec{r}_j\right)^2}{2s}}\,,$$ the average polar order parameter for each cluster size $s$ $$\label{eq:PolarOrderCluster}
P(s)=\left\vert\frac{\sum_{i=1}^s \hat{e}_i}{s} \right\vert\,$$ the tensor order parameter (for 