From the documentation ”Formulation for Observed and Computed Values of Deep Space Network Data Types for Navigation ” expression 4-61 on page 4-42 it can be seen that JPL uses the following expression to account for the effects of relativity under Schwarzschild conditions:
$\frac{d\bar{v}}{dt}=-\frac{GM}{r^2}(1-\frac{4GM}{rc^2}+\frac{v^2}{c^2})\hat{r}+\frac{4GM}{r^2}(\hat{r}\cdot\hat{v})\frac{v^2}{c^2}\hat{r}$
When I use this expression in an integrator it replicates the ”anomalous precession of perihelion” correctly. However, in the Schwarzschild solution in Schwarzschild coordinates you should get the same orbital velocity in a circular orbit as classically. Also, the initial acceleration as you drop an object from rest should be the same as classically (I believe). This JPL expression fails in accomplishing that†. Someone told me that JPL uses isotropic coordinates instead of Schwarzschild coordinates and that this could be an effect of that, but that seems strange to me.
If you use the ”relativistic mass” concept, that works quite well to calculate the relativistic acceleration of a charged particle under influence of the Lorentz force, on gravity you end up with:
$\frac{d\bar{v}}{dt}=-\frac{GM}{r^2}(\hat{r}-\frac{v^2}{c^2}(\hat{r}\cdot\hat{v})\hat{v}) $
This can only generate one third of the perihelion shift, but the expression is better than the JPL expression in the sense that it reproduces correct values for the orbital velocity and the initial acceleration of an object at rest†. By cheating and inserting a factor of three:
$\frac{d\bar{v}}{dt}=-\frac{GM}{r^2}(\hat{r}-3\frac{v^2}{c^2}(\hat{r}\cdot\hat{v})\hat{v}) $
you get an expression that reproduces the correct perihelion shift but also the correct orbital velocity of an object in circular orbit and the initial acceleration of an object at rest.
†The condition for circular motion is $\mathbf{v} \cdot \mathbf{r}=0$, there is no radial part of the motion. Then you set the acceleration terms that do not vanish for $\mathbf{v} \cdot \mathbf{r}=0$ equal to $v^2/r$, the centrifugal acceleration and solve. You see that in the case of no motion,$v=0$, and the case of no radial motion the second and the third expression above reduces to the classical Newtonian gravitational acceleration, which is expected also from the Schwarzschild solution in Schwarzschild coordinates, but the "JPL expression" do not. I would be very happy if someone from JPL could tell me why you are using the first expression above. There is a rudimentary derivation of the expression in the documentation but it is rather high level and not so easy to comprehend.
Note that according to JPL $v = \sqrt{GM/R}$ no longer holds true for a circular orbit but instead you have†:
$v=\sqrt{\frac{GM}{r}\frac{(1-4GM/(rc^2))}{(1-GM/(rc^2))}}$
Also when dropping an object from rest, according to JPL, the acceleration goes as:
$\frac{d\bar{v}}{dt}=-\frac{GM}{r^2}(1-\frac{4GM}{rc^2})\hat{r}$
From this last expression we actually see that JPL, in all of their ephemeris calculations, actually use a small "negative inverse r cube" gravitational term, which is a bit odd.
Questions:
1.Why is JPL using the first expression above and not something similar to the third?
2.What is the correct expression for the orbital velocity for a body in circular motion according to JPL?
3.What is the correct initial acceleration of an object at rest according to JPL?
I would be very happy to get some answers.
Strong field orbits
I did spend a lot of time, see old messy paper, trying to come up with some physical explanation for why, at least in the weak field limit, the third expression above should hold true by experimenting with a "general relativistic relativistic mass" of the type $\gamma(r,v)$ instead of just $\gamma(v)$ but I did not quite succeed. If you insert $\gamma=\frac{1}{\sqrt{1-\frac{2GM}{rc^2}-\frac{v^2}{c^2}}}\frac{1}{\sqrt{1-\frac{2GM}{rc^2}}}$ into $\frac{d(m\gamma\bar{v})}{dt}=-\frac{GMm\gamma}{r^2}\hat{r}$ you end up with $\frac{d\bar{v}}{dt}=-\frac{GM}{r^2}\left(\hat{r}-3\frac{v^2(\hat{r}\cdot\hat{v})\hat{v}}{c^2(1-\frac{2GM}{rc^2})} +\frac{v^4(\hat{r}\cdot\hat{v})\hat{v}}{c^4(1-\frac{2GM}{rc^2})^2}\right)$.
In the strong field limits this expression results in orbits as shown below where the green circle represents the Schwarzschild radius and the red circle represents the radius of the "innermost stable circular orbit" located at a distance of three Schwarzschild radiuses. The result is similar to what is expected from GR.
If you use the JPL-formula in the strong field limit you can get very strange "bouncing" effects as shown below, this is because of the repulsive inverse r-cube term:
This is not at all what is expected from GR. I realized there is a higher order version of the JPL formula that includes an attractive inverse $r^4$ term as well as a repulsive inverse $r^5$ term. Still I think it is very strange to simulate GR by using a repulsive $r^3$ term and I do not really know the reason for why it is common practice to do just that.