Sampling methods to describe photon statistics in large ensembles of quantum emitters

Photon-mediated interactions between quantum emitters can give rise to correlated and directional radiative emission. Many quantum optical phenomena can be characterized by the two-time quantum correlation function g(2)(t, τ), which describes the photon statistics of emitted radiation. However, the critical task of determining g(2)(t, τ) becomes intractable for large emitter ensembles due to the exponential scaling of the Hilbert space dimension with the number of emitters. Here, we analyze and benchmark two approximate numerical sampling methods applicable to emitter arrays embedded within electromagnetic environments, which generally provide upper and lower bounds for g(2)(t, 0). We also introduce corrections to these methods (termed offset corrections) that significantly improve the quality of the predictions. The optimal choice of method depends on the total number of emitters, such that, taken together, the two approaches provide accurate descriptions across a broad range of important regimes. This work, therefore, provides new theoretical tools for studying photon statistics and correlation functions in large ensembles of quantum emitters.