general relativity - What is the exact gravitational force between two masses including relativistic effects?

I was wondering if there is a closed-form formula for the force between two masses $m_1$ and $m_2$ if relativistic effects are included. My understanding is that the classic formula $G \frac{m_1 m_2}{r^2}$ is just an approximation (which is good enough for probably even going to the Moon), but what would the "correct" formula be according to the general theory of relativity? Does a closed formula even exist? For example for an idealized situation of just two spherical masses with homogenous mass distribution?

Answer

Worse than electrodynamics, general relativity is non-linear, in the sense that the field from multiple sources is not just the sum of fields from each isolated source. Even the simple case you are asking about, which is the two-body problem in general relativity, has not been solved exactly.

An even simpler case is the limit $m_2 \to 0$. In that case only $m_1$ affects the geometry of spacetime, and $m_2$ follows a geodesic in that spacetime. This has been solved exactly. The linked article gives details.

[Addendum] To directly answer the original question for the limit $m_2 \to 0$:

Of course if $m_2=0$, then the force between the two masses is $0$. But the thing about classical gravity is that the acceleration due to gravity is independent of mass. (This is the equivalence principle and is actually one of the starting points of the theory of general relativity.) So it still makes sense to ask what would be the acceleration of a negligible mass (aka, a test body) due to the gravity of another mass. The classical answer is $a_2 = \frac{Gm_1}{r^2}$.

In general relativity, the acceleration of a test body due to the gravity of a single spherical, homogenous, non-rotating mass is given exactly by the Schwarzschild solution, which link you can consult for details. The result of which is that

$$\ddot{r} = -\frac{Gm_1}{r^2} + r\dot{\theta}^2 - \frac{3Gm_1}{c^2}\dot{\theta}^2$$ $$\ddot{\theta} = -\frac{2}{r}\dot{r}\dot{\theta}$$ [CORRECTED AND SIMPLIFIED Jan 2]

where $r$ and $\theta$ are polar co-ordinates centred on the gravitating mass, the dots represent differentiation by the proper time of the test body.

So the first term is just the classical radial acceleration $-\frac{Gm_1}{r^2}$. The terms $r\dot{\theta}^2$ and $-\frac{2}{r}\dot{r}\dot{\theta}$ are the classical centrifugal and Coriolis acceleration for polar co-ordinates.

What's not classical is the extra term $\frac{3Gm_1}{c^2}\dot{\theta}^2$. Finally, there is the fact that differentiation is with respect to proper time of the test body. Different test bodies will experience time differently. They can be related by:

$$(1 - \frac{r_s}{r})\dot{t}^2 - \frac{\dot{r}^2}{(1 - \frac{r_s}{r}) c^2} - \frac{r^2\dot{\theta}^2}{c^2} = 1$$

The constant $r_s = \frac{2Gm_1}{c^2}$ is introduced for simplicity. It is called the Schwarzschild radius of the gravitating mass.

Here the co-ordinate $t$ is introduced as a reference time, so $\dot{t}$ is the rate of change of reference time with respect to proper time of the test body. For a distant ($r \to \infty$), stationary ($\dot{r}=0, \dot{\theta}=0$) test body, this becomes $\dot{t} = 1$, so the reference time can be interpreted as the time measured on a distant, stationary clock.

In the classical case of course every body experiences the same time, but it can also be compared to the special relativistic case, where the equation would be:

$$\dot{t}^2 - \frac{\dot{r}^2}{c^2} - \frac{r^2\dot{\theta}^2}{c^2} = 1$$

So what's new in general relativity is the factor $(1 - \frac{r_s}{r})$. To get some idea of scale, for the Earth, $\frac{r_s}{r}$ is about one and a half parts per billion on the surface of the Earth. (Note the Schwarzschild solution is only applicable outside the gravitating body.)

This is exact only for $m_2 \to 0$, but it remains a very good approximation as long as $m_2$ is much smaller than $m_1$, such as for a planet orbiting a star.