  We derive equations for the strongly coupled system of light and dense atomic
ensembles. The formalism includes an arbitrary internal level structure for the
atoms and is not restricted to weak excitation of atoms by light. In the low
light intensity limit for atoms with a single electronic ground state, the full
quantum field-theoretical representation of the model can be solved exactly by
means of classical stochastic electrodynamics simulations for stationary atoms
that represent cold atomic ensembles. Simulations for the optical response of
atoms in a quantum degenerate regime require one to synthesize a stochastic
ensemble of atomic positions that generates the corresponding quantum
statistical position correlations between the atoms. In the case of multiple
ground levels or at light intensities where saturation becomes important, the
classical simulations require approximations that neglect quantum fluctuations
between the levels. We show how the model is extended to incorporate
corrections due to quantum fluctuations that result from virtual scattering
processes. In the low light intensity limit we illustrate the simulations in a
system of atoms in a Mott-insulator state in a 2D optical lattice, where
recurrent scattering of light induces strong interatomic correlations. These
correlations result in collective many-atom subradiant and superradiant states
and a strong dependence of the response on the spatial confinement within the
lattice sites.
