Is there a rigorous proof regarding the non-linear stability of the $L_4$ and $L_5$ Lagrange points?

I have found that many proofs regarding the stability of the $$L_4$$ and $$L_5$$ Lagrange points are based on linear approximations of the equations of motion near these points. However, from a dynamical systems perspective, Lagrange points are essentially fixed points of the phase flow. Linear approximations at fixed points yield a linear system of equations, such as $$\dot{\phi}=A\phi$$, where $$\phi$$ is a coordinate in phase space and $$A$$ is a constant matrix. If all eigenvalues of $$A$$ have negative real parts, then the fixed point is stable. If $$A$$ has an eigenvalue with a positive real part, then the fixed point is unstable, which corresponds to the case of $$L_1$$ to $$L_3$$ points. However, if all eigenvalues of $$A$$ are purely imaginary, we cannot determine the stability of the fixed point from the linear approximation of the phase flow. This is precisely the case for $$L_4$$ and $$L_5$$ points.

Many "proofs" I have come across only argue that the eigenvalues are imaginary, without considering higher-order approximations of the phase flow, cf. e.g. this Phys.SE post. I would like to inquire if such treatment is rigorous, and if there exists a rigorous proof regarding the stability of the $$L_4$$ and $$L_5$$ Lagrange points.

To illustrate the details, here is an example from "An Introduction to Dynamical Systems: Discrete and Continuous" by R. C. Robinson:

Consider a 2-D system: \begin{aligned} \dot{x} &= -y+ax(x^2+y^2)\\ \dot{y} &= x+ay(x^2+y^2) \end{aligned} The origin is a fixed point of this system, and through linear approximation, we can determine that $$\begin{bmatrix} \delta\dot x\\ \delta\dot y \end{bmatrix}=\begin{bmatrix} 0 & -1\\ 1 & 0 \end{bmatrix}\begin{bmatrix} \delta x\\ \delta y \end{bmatrix}$$ The eigenvalues of this linear approximation matrix are $$\pm\mathrm{i}$$, which are purely imaginary. However, if we consider this system in polar coordinates, using $$r^2=x^2+y^2$$, we have: \begin{aligned} &r\dot r = x\dot x+y\dot y=x(-y+axr^2)+y(x+ayr^2)=ar^4\\ \Rightarrow\quad& \dot r= ar^3 \end{aligned} It can be seen that when $$a>0$$, the origin is a repelling point and unstable; only when $$a\le0$$, the origin is a stable point.

• The "approximations" are completely rigorous from the point of view of establishing that stable orbits exist around the Lagrange points. The question isn't if you can orbit the Lagrange points arbitrarily far from the Lagrange points - of course you can't. The question is whether or not there exist stable orbits very close to the Lagrange points. The linear approximation becomes better and better as you consider orbits closer and closer to the Lagrange points. If you aren't satisfied with the accuracy of the approximation, then consider orbits even closer to the Lagrange points. Oct 17, 2023 at 9:13

I found part of the answer to this question in a doctoral thesis: stability of Lagrangian Points. This thesis claims that the stability of the $$L_4$$ and $$L_5$$ points in the planar case can be rigorously proven using the KAM theory. However, for the three-dimensional case, even with the use of KAM theory, the stability of the $$L_4$$ and $$L_5$$ points cannot be proven. More importantly, the thesis mentions in the conclusion: "In addition, the infinitesimal particle is thought to leave from a small region near these equilibria, through a phenomenon currently known as Arnold's diffusion. Nevertheless, the time that it may take for this particle to leave from an equilibrium point might be of the order of the age of the solar system." If this statement is correct, then mathematically speaking, the $$L_4$$ and $$L_5$$ points in three-dimensional space are unstable (in the sense of Lyapunov stability).