  Coalescing massive black hole binaries, formed during galaxy mergers, are
expected to be a primary source of low frequency gravitational waves. Yet in
isolated gas-free spherical stellar systems, the hardening of the binary stalls
at parsec-scale separations owing to the inefficiency of relaxation-driven
loss-cone refilling. Repopulation via collisionless orbit diffusion in triaxial
systems is more efficient, but published simulation results are contradictory.
While sustained hardening has been reported in simulations of galaxy mergers
with $N \sim 10^6$ stars and in early simulations of rotating models, in
isolated non-rotating triaxial models the hardening rate continues to fall with
increasing N, a signature of spurious two-body relaxation. We present a novel
approach for studying loss cone repopulation in galactic nuclei. Since loss
cone repopulation in triaxial systems owes to orbit diffusion, it is a purely
collisionless phenomenon and can be studied with an approximated force
calculation technique, provided the force errors are well behaved and
sufficiently small. We achieve this using an accurate fast multipole method and
define a proxy for the hardening rate that depends only on stellar angular
momenta. We find that the loss cone is efficiently replenished even in very
mildly triaxial models (with axis ratios 1 : 0.9 : 0.8). Such triaxiality is
unavoidable following galactic mergers and can drive binaries into the
gravitational wave regime. We conclude that there is no final parsec problem.
