How is relativistic $N$-body simulation possible without knowing the entire histories of each particle? In $N$-body simulation you need to know the positions of the particles in order to calculate the force between them. The new velocity of each particle can then be calculated given a simulation timestep $dt$.
If the gravitational interaction propagates at the speed of light, do we need to specify the force between particles given their retarded positions (not instantaneous positions)? If not, why not?
Note that I am interested in special relativity only (not GR), I am more concerned with what happens when bodies move quickly or are separated by distances larger than the distance light travels in one timestep.
I am also interested in whether or not the electric force between charges is any different in this respect.
 A: Steve Carlip has answered in a more intuitive way the question of where the EM field and where the gravitational field "points" when the source is moving.  He a short version on John Baez's physica FAQ
http://math.ucr.edu/home/baez/physics/Relativity/GR/grav_speed.html
and a more technical paper https://arxiv.org/ftp/arxiv/papers/1912/1912.05439.pdf
Essentially the short answer is EM radiation's lowest mode is dipole so it requires a second derivative of the motion to radiate.  So the field at a distance, points at the position extrapolated by assuming constant velocity over the propagation time.  Gravity radiation is quadrapole and so requires a third derivative of the motion to radiate.  Then the field points at a position extrapolated by assuming constant acceleration over the propagation time to the distant point.
A: I believe that you are not being coherent.
As far as I know, you can only choose to treat gravity on either a non-relativistic regime, where you totally neglects any relativistic effect, and thus you have newtonian gravity, or you have to use the full machinery of GR. I believe that there is no 'Special-Relativistic gravity' or 'Semi-Relativistic gravity'. In the case of newtonian gravity you would specify on the instantaneous positions, but because of the regime that you are on, that is not a problem in itself.
If you particles can even have the chance to move on a time-step faster than light, you have long departed the regime of validity of newtonian gravity, so, it's not the implementation that is wrong, but the modeling of the problem.
In the case of EM you have that Maxwell's equation are covariant with respect to lorentz transformations, so there should be no problem if you use covariant equations for the source $(\rho,\vec j)$. 
@user34168
I don't know. That is a very good question, but my first thought is that since gravity couples to the full Stress-Energy Tensor, and the idea to get newtonian gravity is to do $T^00 \approx \rho c^2 $ and the rest you ignore, this surelly fails when you have relativistic speeds, even though the metric could be safelly aproximated.
I believe that you should be able to to something like retarded potentials + gravitational radiation corrections, but I really didn't do the calculation so I'm not going to go on speculating.
About the Coulumb x Liénard–Wiechert potential. Yes, this is the first idea, but as far as I remember what happens is that you have a bunch of point charges, you calculate their LW potential, then you apply that potential on something else, but you don't calculate the back reation of the target on the source. Also, I don't know if you can readly extend the LW potentials to non-point charges, but I believe not.
Edit:
 Corrected ambiguity
Edit 2:
 @Ruslan, thanks for noting the typo
A: You have essentially discovered the need for field theories.  Instead of knowing arbitrarily long sections of the particle's past world lines, the alternative is to specify the field variables on a time-slice.
