Chapter 4 Problem 4.11: The precession of the Perihelion of Mercury

计算物理

---

Name:陈锋
Student Number: 2013301020145
May 5, 2016

Abstract

The author use fourth-order Runge-Kutta method to solve this problem and the algorithm is contained in ode.py

Introduction

Excepting the influence of Jupiter and Saturn, there is an additive unexplained 43 arcseconds of precession of Mercury. Using general relativity, it can be explained theoretically. However, the precession rate is fairly small, so we have to design our simulation with that in mind.
The force law predicted by general relativity is:

$$\begin{equation} F_G \approx -\frac{GM_SM_M}{r^2}(1+\frac{\alpha}{r^2}) \end{equation}$$
where $M_M$ is the mass of Mercury and $\alpha \approx 1.1 \times 10^{-8} AU^2$. In order to observe the obvious precession of Mercury, we are willing to increase $\alpha$ to 0.01.

However, if the coefficient $\alpha$ is too large, the trajectory of Mercury will deviate from ellipse, which will bring troubles to us.

The Equation of motion

The orbital trajectory fro a body of reduced mass $\mu$ is given in polar coordinates by

$$\begin{equation} \frac{d^2}{d\theta^2}(\frac{1}{r}) + \frac{1}{r} = -\frac{\mu r^2}{L^2}F(r) \end{equation}$$
where $L = \mu r^2 \dot{\theta}$ is the angular momentum and $F(r)$ is the force acting on the body.
According to general relativity:
$$\begin{equation} \frac{d^2}{d\theta^2}(\frac{1}{r}) + \frac{1}{r} = \frac{\mu GM_SM_M}{L^2}(1+\frac{\alpha}{r^2}) \end{equation}$$
Since $\alpha$ is very small, the solution of equation (3) must be similar to that of equation (2) if $F(r) = - GM_SM_M/r^2$. it means that:
$$\begin{equation} r_{max} = a(1+e) \\ r_{min} = a(1-e) \\ L = \sqrt{\mu G M_S M_M a (1 - e^2)} \end{equation}$$ where e is the eccentricity and a is the major semi-axis.
Thus, the equation becomes:
$$\begin{equation} \frac{d^2}{d\theta^2}(\frac{1}{r}) + \frac{1}{r} = \frac{1}{a(1-e^2)}(1+\frac{\alpha}{r^2}) \end{equation}$$
Transferring it into differential equations:
$$\begin{equation} \frac{dr}{d\theta} = v \\ \frac{dv}{d\theta} = \frac{2}{r}v^2 + r - K(1+\frac{\alpha}{r^2})r^2 \end{equation}$$
where $K = 1/a(1-e^2)$. Via these euqations, we can use numerical method to solve this problem.

The Accuracy of Numerical Method

Firstly, if we doesn't consider the effect of general relativity, the trajectory can be solve analytically. Thus, in this situation, we can check the precision of our numerical method.

Figure 1: The solid line is the result of numerical method and the dash line is that of analytical method. Since I use fourth-order Runge-Kutta method, the numerical results are nearly the same as the analytical one, which is what we expect. The unit of r is AU.

The Precession of diverse elliptical orbits with different eccentricities

The simulation of percession

Figure 2: Simulated orbit for Mercury orbiting the sun. The force law (1) was used, with $\alpha = 0.01$. The programme was stopped after several orbits.

Fixing the perihelion

Since we know that:

$$\begin{equation} K = \frac{1}{a(1-e^2)} = \frac{1}{r_{min}(1+e)} \\ r_{max} = \frac{1+e}{1+e}r_{min} \end{equation}$$
If the influence of general relativity is not so large, we can regard the orbit of Mercury as a ellipse. Then, using the equation (7), we can fix the perihelion when changing the value of eccentricity.

The influence on the shapes of the orbits

Setting $\alpha = 0.001$ and $r_{min} = 0.3 AU$

Figure 3: As the increase of eccentricity, the furthest distance between Mercury and the sun becomes larger.

The influence on angular velocity

Figure 4: As the increase of eccentricity, the angular velocity of percession decrease. It is intuitive. Because when the eccentricity is very large, the oribit will become so long that it takes a very long period to circulate.

Fitting result:

Linear model	f(x) = a*exp(b*x)
Coefficient: a (with 95% confidence bounds)	21.13 (19.15, 23.11)
Coefficient: b (with 95% confidence bounds)	-3.67 (-4.076, -3.264)

Goodness of fit:

SSE	1.415
R-square	0.9925
Adjusted R-square	0.9914
RMSE	0.4496

Conclusion

As the eccentricity increase, the orbit of Mercury becomes longer and the angular velocity of percession dramatically decrease.

Acknowledgement

感谢王世兴同学指出计算不同离心率对应的进动角速度时作者所犯的错误。

Reference

1.Giodano, N.J., Nakanishi, H. Computational Physics. Tsinghua University Press, December 2007.

Programme Code

• Chapter 4 Problem 4.11
• ode.py