While [Eq. \[eq:FTMS\_output\_field\]]{} and [Eq. \[eq:FTMS\_coupled\_dipoles\]]{} are formally exact, solving a number of equations that explicitly scales with the number of atoms and that depends on the details of atomic positions is not a particularly convenient way to calculate the index or other optical properties. Historically, this fostered the development of simplified theories for the macroscopic response, such as the Drude-Lorentz model [@Jackson1998ClassicalElectrodynamics] or equivalently the Maxwell-Bloch (MB) equations [@Grynberg2010IntroductionOptics], where the discreteness of atoms is replaced by a smooth medium of density $N/V$ ([Fig. \[fig:scatt\_cross\_sec\]]{}-c). The resulting index depends on the product of density and single-atom polarizability, [ $$\begin{array}{c}

\label{eq:FTMS_index_prediction_MB}
n_{\text{MB}}(\Delta)=\sqrt{1+\dfrac{N}{V}\alpha_0(\omega)}= \sqrt{1+ \dfrac{3\pi\eta}{-\Delta - i/2}},

\end{array}$$ ]{} where we defined the dimensionless density $\eta \equiv N/(Vk_0^3)$. Notably, for an optimum detuning, the maximum real part of the index scales as $\sim\sqrt{\eta}$.

While the MB equations ignore multiple scattering, the Lorentz-Lorenz (LL) or the equivalent Clausius-Mossotti model is one well-known approach to approximate its effects, still within the smooth density approximation. Given any atom located at ${\mathbf{r}}_0$, the model approximates the neighbouring atoms as a smooth dielectric medium with a small spherical exclusion around ${\mathbf{r}}_0$ [@Jackson1998ClassicalElectrodynamics]. The resulting local field correction gives an index that satisfies the equation $(n_{\text{LL}}^2-1)/(n_{\text{LL}}^2+2) = 