  A particle method for reproducing the phase space of collisionless stellar
systems is described. The key idea originates in Liouville's theorem which
states that the distribution function (DF) at time t can be derived from
tracing necessary orbits back to t=0. To make this procedure feasible, a
self-consistent field (SCF) method for solving Poisson's equation is adopted to
compute the orbits of arbitrary stars. As an example, for the violent
relaxation of a uniform-density sphere, the phase-space evolution which the
current method generates is compared to that obtained with a phase-space method
for integrating the collisionless Boltzmann equation, on the assumption of
spherical symmetry. Then, excellent agreement is found between the two methods
if an optimal basis set for the SCF technique is chosen. Since this
reproduction method requires only the functional form of initial DFs but needs
no assumptions about symmetry of the system, the success in reproducing the
phase-space evolution implies that there would be no need of directly solving
the collisionless Boltzmann equation in order to access phase space even for
systems without any special symmetries. The effects of basis sets used in SCF
simulations on the reproduced phase space are also discussed.
