frequencies $\omega_{i,j}$ with the new values $\omega_{\pm}$ (second panel of [Fig. \[fig:RG\_plots\]]{}-b). While the resulting modes are in principle delocalized between atoms $i,j$, to facilitate the RG, we randomly assign $\omega_{+}$ to either atom $i$ or $j$, while $\omega_{-}$ is then assigned to the other atom (see SI on the issue of replacing atoms $i,j$ with two new atoms placed at the midpoint of the original locations). This new system is described by a renormalized interaction matrix $\mathcal M_{\text{eff}}=\text{diag}(\boldsymbol{\omega}_{\text{eff}})-\tilde G_{\text{eff}}$, where $\boldsymbol{\omega}_{\text{eff}}=(\omega_1,\dots,\omega_{+}, \dots, \omega_{-},\dots \omega_N)$ contains the two renormalized resonance frequencies, and where $\tilde G_{\text{eff}}$ includes the new set of allowed near-field interactions $\mathcal L_{\text{eff}}$, which both forbid the renormalized pair from interacting again (i.e. $\mathcal L^{\text{eff}}_{ij}=0$) and prevent any backflow of the RG process (see Appendix A for more details). The RG process can be iteratively repeated by identifying, at each step, the most strongly interacting pairs, and ends once $\mathcal K_{ij} \leq \mathcal K_{\text{cut-off}}\sim 1$, i.e. when all strong near-field interactions have been removed. In the numerics presented here, we take a cutoff parameter of $\mathcal K_{\text{cut-off}}= 1$. Other choices result in minor quantitative corrections, while the overall conclusions remain the same. The final result, as suggested in 