How do simulations of galaxies or star systems account for the propagation delay of gravity?

Question

I've made simple $N$-body gravity simulations as a programming and physics exercise in the past. But they've never accounted for speed-of-light time delay.

The way I'm thinking about this is illustrated by an example: If I have two stars, A and B, that orbit each other at say 10 light minutes distance, but my time-step in the simulation is 1 minute. So for star A, when computing the force acting on it from B, it should use the position of B from 10 steps ago (10 light minutes, 1 minute time step) rather than it's current position. And the reverse for star B's perspective.

Is this even the right way to think about this? Or is there a modification of the numerical solution to the gravity equations that would be more accurate? Would have to code all previous positions for each entity?

I've also had the thought that if the timestep is on the order of D/c where D is total size of the simulation then the timestep is the same as the time delay and I guess the problem would just wash out.

Answer 1

I am evading your question here and instead answer by saying: You only need to consider this if over the course of the propagation of gravity your positions change notably. Typically in galaxies and stars that is not the case, since speeds are very non-relativistic, and the distances are large, i.e. the change in the force vector is tiny and the Newtonian approximation holds. Then, going beyond that, for the solar system, you would go to the so-called Post-Newtonian Approximation, which is the proper series expansion of General Relativity in the appropriate limit. If you only take a finite speed of gravity into account but otherwise still calculate with a Newtonian force then you are likely to miss out on more important effects.