<p>In Newtonian gravity, a self-gravitating collisionless gas around a massive object such as a star or a planet is modeled via the Vlasov–Poisson system with an external Kepler potential. The presence of this attractive potential allows for bounded trajectories along which the gas neither falls in towards the object nor escape to infinity. We study this system focusing on the regime with bounded trajectories. First, we prove quantitative linear phase mixing estimates in three dimensions outside symmetry. Second, our main result is a long-time nonlinear phase mixing theorem for spherically symmetric data with finite regularity. The mechanism is phenomenologically similar to Landau damping on a torus and our result applies to the same time scale (modulo logarithms) as the known results on Landau damping with finite regularity. However, in contrast with Landau damping, we need to contend with weaker linear estimates as well as use a system of dynamically defined action angle variables.</p>