  Dark matter density is formally infinite at the location of caustic surfaces,
where dark matter sheet folds in phase space. The caustics separate
multi-stream regions with different number of streams. Volume elements change
the parity by turning inside out when passing through the caustic stage. Being
measure-zero structures, identification of caustics via matter density fields
is usually restricted to fine-grained simulations. Instead a generic purely
geometric algorithm can be employed to identify caustics directly by using
triangulation of Lagrangian sub-manifold x(q, t) where x and q are Eulerian and
Lagrangian coordinates obtained in N-body simulations. The caustic surfaces are
approximated by a set of triangles with vertices being the particles in the
simulation. It is demonstrated that finding a dark matter halo is quite
feasible by building its owtermost convex caustic. Neither more specific
assumptions about the geometry of the boundary nor ad hoc parameters are
needed. The halo boundary in our idealized but undoubtably generic simulation
is neither spherical nor ellipsoidal but rather remarkably asymmetrical. The
analysis of the kinetic and potential energies of individual particles and the
halo as a whole along with an examination of the two-dimensional phase space
has shown that the halo is gravitationally bound.
