  We use the Delaunay Tessellation Field Estimator (DTFE) to study the
one-point density distribution functions of the Millennium (MS) and
Millennium-II (MS-II) simulations. The DTFE technique is based directly on the
particle positions, without requiring any type of smoothing or analysis grid,
thereby providing high sensitivity to all non-linear structures resolved by the
simulations. In order to identify the detailed origin of the shape of the
one-point density probability distribution function (PDF), we decompose the
simulation particles according to the mass of their host FoF halos, and examine
the contributions of different halo mass ranges to the global density PDF. We
model the one-point distribution of the FoF halos in each halo mass bin with a
set of Monte Carlo realizations of idealized NFW dark matter halos, finding
that this reproduces the measurements from the N-body simulations reasonably
well, except for a small excess present in simulation results. This excess
increases with increasing halo mass. We show that its origin lies in
substructure, which becomes progressively more abundant and better resolved in
more massive dark matter halos. We demonstrate that the high density tail of
the one-point distribution function in less massive halos is severely affected
by the gravitational softening length and the mass resolution. In particular,
we find these two parameters to be more important for an accurate measurement
of the density PDF than the simulated volume. Combining our results from
individual halo mass bins we find that the part of the one-point density PDF
originating from collapsed halos can nevertheless be quite well described by a
simple superposition of a set of NFW halos with the expected cosmological
abundance over the resolved mass range. The transition region to the
low-density unbound material is however not well captured by such an analytic
halo model.
