at the heart of the RG approach for the many-atom case.

We now discuss how strong, coherent $1/\rho_{ij}^3$ near-field interactions in a many-atom system can be treated, by successively replacing strongly interacting pairs by optically equivalent, non-interacting atoms. Here, we will focus on the main conceptual steps of our RG scheme, while full details can be found in the Appendices. We anticipate that the scheme involves atoms with different renormalized resonant frequencies $\omega_i$, which can either interact, or not, through the near-field coupling, depending on the previous RG steps. The normal modes of such a system are given by the eigenstates of the generalized $N \times N$ matrix $\mathcal M= \text{diag}( \boldsymbol \omega )-\tilde G$, where the elements $\tilde G_{ij}$ are defined as $\tilde G_{ij}= \mathcal L_{ij} G_{ij}^{\text{near}} + (G_{ij}-G_{ij}^{\text{near}})$. Here, $\text{diag}( \boldsymbol \omega )$ is a diagonal matrix containing the individual resonance frequencies $\boldsymbol \omega =(\omega_1,\dots,\omega_N)$, while $\mathcal L_{ij}=1$ or $0$ dictates whether pair $i,j$ is allowed to interact via the near field ($\mathcal L_{ij}=1$ for all pairs at the beginning of the RG process). In three dimensions, the $1/\rho^3$ scaling of the near-field interaction implies that if an atom has a particularly close-by nearest neighbor, this pair will 