  Time-symmetric integration schemes share with symplectic schemes the property
that their energy errors show a much better behavior than is the case for
generic integration schemes. Allowing adaptive time steps typically leads to a
loss of symplecticity. In contrast, time symmetry can be easily maintained, at
least for a continuous choice of time step size. In large-scale N-body
simulations, however, one often uses block time steps, where all time steps are
forced to take on values as powers of two. This greatly facilitates
parallelization, and hence code efficiency. Straightforward implementation of
time-symmetry, translated to block time steps, faces significant hurdles. For
example, iteration can lead to oscillatory behavior, and even when such
behavior is suppressed, energy errors show a linear drift in time. We present
an approach that circumvents these problems.
