  We present a new and completely general technique for calculating the
fine-grained phase-space structure of dark matter throughout the Galactic halo.
Our goal is to understand this structure on the scales relevant for direct and
indirect detection experiments. Our method is based on evaluating the geodesic
deviation equation along the trajectories of individual DM particles. It
requires no assumptions about the symmetry or stationarity of the halo
formation process. In this paper we study general static potentials which
exhibit more complex behaviour than the separable potentials studied
previously. For ellipsoidal logarithmic potentials with a core, phase mixing is
sensitive to the resonance structure, as indicated by the number of independent
orbital frequencies. Regions of chaotic mixing can be identified by the very
rapid decrease in the real space density of the associated dark matter streams.
We also study the evolution of stream density in ellipsoidal NFW halos with
radially varying isopotential shape, showing that if such a model is applied to
the Galactic halo, at least $10^5$ streams are expected near the Sun. The most
novel aspect of our approach is that general non-static systems can be studied
through implementation in a cosmological N-body code. Such an implementation
allows a robust and accurate evaluation of the enhancements in annihilation
radiation due to fine-scale structure such as caustics. We embed the scheme in
the current state-of-the-art code GADGET-3 and present tests which demonstrate
that N-body discreteness effects can be kept under control in realistic
configurations.
