  Modeling self-gravity of collisionless fluids (e.g. ensembles of dark matter,
stars, black holes, dust, planetary bodies) in simulations is challenging and
requires some force softening. It is often desirable to allow softenings to
evolve adaptively, in any high-dynamic range simulation, but this poses unique
challenges of consistency, conservation, and accuracy, especially in
multi-physics simulations where species with different softening laws may
interact. We therefore derive a generalized form of the energy-and-momentum
conserving gravitational equations of motion, applicable to arbitrary rules
used to determine the force softening, together with consistent associated
timestep criteria, interaction terms between species with different softening
laws, and arbitrary maximum/minimum softenings. We also derive new methods to
maintain better accuracy and conservation when symmetrizing forces between
particles. We review and extend previously-discussed adaptive softening schemes
based on the local neighbor particle density, and present several new schemes
for scaling the softening with properties of the gravitational field, i.e. the
potential or acceleration or tidal tensor. We show that the tidal softening
scheme not only represents a physically-motivated, translation and Galilean
invariant and equivalence-principle respecting (and therefore conservative)
method, but imposes negligible timestep or other computational penalties,
ensures that pairwise two-body scattering is small compared to smooth
background forces, and can resolve outstanding challenges in properly capturing
tidal disruption of substructures (minimizing artificial destruction) while
also avoiding excessive N-body heating. We make all of this public in the GIZMO
code.
