  $N$-body simulations study the dynamics of $N$ particles under the influence
of mutual long-distant forces such as gravity. In practice, $N$-body codes will
violate Newton's third law if they use either an approximate Poisson solver or
individual timesteps. In this study, we construct a novel $N$-body scheme by
combining a fast multipole method (FMM) based Poisson solver and a time
integrator using a hierarchical Hamiltonian splitting (HHS) technique. We test
our implementation for collision-less systems using several problems in
galactic dynamics. As a result of the momentum conserving nature of these two
key components, the new $N$-body scheme is also momentum conserving. Moreover,
we can fully utilize the $\mathcal O(\textit N)$ complexity of FMM with the
integrator. With the restored force symmetry, we can improve both angular
momentum conservation and energy conservation substantially. The new scheme
will be suitable for many applications in galactic dynamics and structure
formation. Our implementation, in the code Taichi, is publicly available at
https://bitbucket.org/qirong_zhu/taichi_public/.
