As mentioned in the introduction, the so-called Runge-Lenz vector was discovered by Hermann, transmitted by Bernoulli, used intensively by Laplace, exhumed by Runge at the beginning of the 20^th century and used by Lenz and Pauli in quantum mechanics in order to study the hydrogen atom. Although the vector should therefore be called Hermann-Bernoulli-Laplace-Runge-Lenz vector [4,5,21], we will in the following, for simplicity and because it is the most common way it is referred to in the literature, retain the denomination “Runge-Lenz” vector. The Runge-Lenz vector reads [22,23]: where $\mathbf{v}$ is the velocity of our system (mass m), $\mathbf{L}=\mu \mathbf{r}\wedge \mathbf{v}$ the angular momentum, and ${\mathbf{u}}_r=\mathbf{r}/r$ the radial unit vector. Throughout the paper, the symbol ∧ denotes the vector product:

The equation of motion is and one has

It is also possible to obtain it more directly. Using the invariance of angular momentum which gives and taking the vector product of the latter equality by $\mathbf{L}$, one has giving directly the invariance of the Runge-Lenz vector. The projection of the Runge-Lenz vector on the big axis of the ellipse gives the equation of the conic. Calculating the scalar product $\mathbf{A}.\mathbf{r}$, one gets where we have used The polar equation of the orbit is then [24,25]: with $p=C^2/k=L^2/\left(\kappa {\mu }^2\right)$, $C=L/\mu$ being the areal constant (twice the areal velocity, Kepler’s second law) and $e=A/\left(k\mu \right)=A/\alpha$ the eccentricity of the ellipse. The energy of the motion can be obtained by calculating the norm of the Runge-Lenz vector Using the definitions of p and e, as well as the relation $p=a(1-e^2)$, we deduce without using the Binet formulas [26] (see Section 2.2). There are other invariant vectors, which are related to the Runge-Lenz vector.

In celestial mechanics, the eccentricity vector of a Kepler orbit [27,28] is the dimensionless vector with direction pointing from apoapsis to periapsis (respectively farthest and nearest points in the orbit of a planetary body about its primary body to the center of force) and with magnitude equal to the orbit’s scalar eccentricity. For Kepler orbits, the eccentricity vector is a constant of motion. It is of great interest in the analysis of almost circular orbits, as perturbing (non-Keplerian) forces on an actual orbit will cause the osculating eccentricity vector to change continuously as opposed to the eccentricity and argument of periapsis parameters for which eccentricity zero (circular orbit) corresponds to a singularity. The eccentricity vector is which follows immediately from the vector identity (10). One has with $\mathbf{p}=\mu \mathbf{v}$. The Runge-Lenz vector $\mathbf{A}$ is directed from the centre of force towards the periapsis, making it easy to visualize (Lenz and Pauli called it the “Achsenvektor”) [29]. However, its physical dimension (SI unit kg m${ }^3$ s${ }^{-2}$) does not endow it with any familiar meaning; the alternative invariant $\mathbf{A}/m$ is hardly better. If $\mathbf{e}$ is used instead of $\mathbf{A}$ in the derivation of the orbit equation, a simple equation appears, with an immediate physical meaning. This is why it is often considered that e should be preferred over A for didactic applications. In the literature, the only “eccentricity-like” vector referred to seems to be the one defined by Hamilton, $\mathbf{L}\wedge \mathbf{A}/L^2$: it is co-linear to $\mathbf{e}$, but has the dimension of a velocity or of a momentum according to the alternative definition of [30]. It does not coincide with a particular value of v or $\mu v$, and thus has no special significance. More importantly, while equations of motion have been obtained with all vectors, each one has its own mathematical simplicity (complexity) and physical meaning. Finally, the eccentricity vector naturally arises in solving the dynamical equation at an elementary level, without resorting to Lagrangian formalism or group theory. The eccentricity vector also only arises in the case of an inverse-square law of force.

The vector product of $\mathbf{L}$ and $\mathbf{A}$ is an invariant. It is thus possible to deduce the invariant Hamilton $\mathbf{S}$ vector The Laplace vector indicates the direction of the periapsis. In the same way, the Hamilton vector is parallel to the velocity of the system at the periapsis [31,32,33,34].

The Bertrand theorem [21,35,36,37,38,39,40,41,42,43,44,45] states that the orbit is closed only if the potential is in $1/r$ or harmonic $r^2$. In the Keplerian case, there are two privileged directions: the one relating the attractive center to the periapsis of the orbit, and the perpendicular direction. The two directions play different roles. This justifies the existence of the vector $\mathbf{A}$. As mentioned above, the vector $\mathbf{A}$ is parallel to the big axis of the ellipse and indicating the direction of the periapsis. In the case of an harmonic ($r^2$) potential, there are also two privileged directions. The Binet formulas are, with $u=1/r$:

where $a_r$ is the acceleration. Let us consider a force proportional to $r^n$ [46]: and a small perturbation $ϵ$ to the circular trajectory of radius $r_0$ (corresponding to $d^{2}u/d{\theta }^2=0$): One has where $ϵ\ll 1$. Using the second Binet formula, one obtains where we have set $2\beta =-(n+2)$ and $6\gamma =(n+2)(n+1)$. We have to solve an equation describing a non-linear oscillation of variable $ϵ$. For that purpose, we can resort to a perturbative expansion of the solution and evaluate the apsidal angle of the orbits in the neighborhood of the circular one. At first order of approximation, in the vicinity of the circular trajectory, one has and thus with $p^2=n+3$. At the next (second) order, we get yielding with $G=(n+2)/4$ and $A=-(n+2)/12$. Finally, at the third order, we look for solutions of the type satisfying the equation if and only if where we have set as well as The solution at third order, around the circular trajectory, is of the form $ϵ=Mcos\left(q\theta \right)+O\left(M^2\right)$, where and the apsidal angle, equal to $\pi /q$ can not be commensurable with $\pi$$\forall M\ll 1$, except if $n=1$ or $n=-2$. Therefore, the orbits of a point subject to a central force proportional to the inverse of the power n of the distance between the active center and the point are closed if and only if $n=1$ or $n=-2$.

The resemblance between the Kepler problem of planetary motion and the charged particle in Coulombic potential motivates us to understand symmetries possessed by such systems. The Hamiltonian can be written formally as where $\mu$ is the reduced mass $Mm/(M+m)$, $\alpha =GMm=G\mu (M+m)=\kappa \mu$ (for Newton potential, see above) and $\alpha =Z{e}^2$ (for Coulomb potential with $4\pi {ϵ}_0$ set to unity). With $\mathbf{r}$, $\mathbf{p}$, and $\mathbf{L}$ as already established quantum mechanical operators, the classical Lenz vector $\mathbf{A}$ can be translated directly into quantum mechanics (vector $\mathbf{M}$). There is only a bit of subtlety involved in its construction. We must take care in defining the cross product (see Equation (3). We notice that $\mathbf{p}\wedge \mathbf{L}\ne -\mathbf{L}\wedge \mathbf{p}$. The classical Runge-Lenz operator is not hermitian. Following the derivation of Ref. [47], we redefine $\mathbf{M}$ as a symmetric average. The quantum-mechanical operator $\mathbf{M}$, extension of the classical Runge-Lenz vector, defined as follows As expected, this vector commutes with the Hamiltonian, $\left[\mathbf{M},H\right]=0$ and defines a conserved quantity in the quantum mechanical analog of the classical Kepler problem. We find that and From here, we look to understand how this new commuting operator relates to the physical system and all the existing conserved quantities. We examine the system in light of well-established geometrical symmetries.

The three components of $\mathbf{M}$, like the three components of $\mathbf{L}$, can be treated as generators of some infinitesimal transformations. With the goal of exploring the algebraic structure of these new generators, we work out their commutation relations (there are thirty-six in all). Some are known; indeed, $\left[L_i,L_j\right]$ are given by for $i,j,k=1,2,3$, where ${ϵ}_{ijk}$ is the Levi-Civita symbol and ${\delta }_{ij}$ the Kronecker symbol. The others are and The first set of these commutators (38) shows that $\mathbf{M}$ is a vector. $M_i$ commutes with $L_i$ since $L_i$ generates rotation about the i axis, which has no effect on a vector pointed along that axis. On the other hand, if j is perpendicular to i, $L_j$ induces a change in the direction of $M_i$. This is characterized by the second commutation relation in Equation (38), which implies that the direction of this change in $M_j$ should be orthogonal to both $M_j$ and $L_i$, in this case, along $-M_k$. The relation in Equation (39) breaks the closed algebra of the L and M operators together, since H is involved. We now make a useful substitution. Considering a degenerate subspace (manifold) of H with energy eigenvalue E, we can replace the Hamiltonian with the resulting energy eigenvalue, and define ${\mathbf{M}}^{\prime }$ such that The commutation relations remain unchanged with M replaced everywhere by $M^{\prime }$, with the exception of Equation (39), which now becomes signature of a closed algebraic system.

We can relabel $L_x$, $L_y$, and $L_z$ into $L_1$, $L_2$, and $L_3$. Substituting 1, 2 and 3 for x, y and z allows for a natural extension beyond three dimensions. The fact that the commutators $\left[r_i,p_j\right]$ are zero only when $i\ne j$ implies that the ${\mathbf{M}}^{\prime }$ matrices are merely generalized angular momentum operators in the fourth dimension. In four dimensions, we have $\left(\begin{array}{c}4\\ 2\end{array}\right)=6$ orthogonal angular momentum generators, where $\left(\begin{array}{c}n\\ p\end{array}\right)=n!/p\S (n-p)!$ is the usual binomial coefficient. Let us make the following associations: which can be summarized as We add a fictitious coordinate $\omega$, so that the rotation in four-dimensional space is given by In order to obtain the energy, we can exploit the $SO\left(4\right)$ symmetry. It is easy to verify that the commutation relations in Section 4.3 still hold. In four dimensions, a degree of decoupling emerges: the maximal commuting set of operators includes two angular momentum operators at the same time: $L_{12}$ with $L_{34}$, etc. Such a decoupling stems from the property that the two operators generate rotations in orthogonal, non-intersecting planes, defined by completely disjoint sets of vectors. This fact alludes to the structure of the group generated by the operators: the rotational group in four dimensions, $SO\left(4\right)$, of rank 2, which contains, as a subgroup, the rank 1 $SO\left(3\right)$ rotational group. We have therefore generalized the closed algebra of $\mathbf{L}$ and $\mathbf{M}$ operators. The remarkable insight here is that a dynamical symmetry (as physically evidenced by eccentricities in closed orbits) is a mere artifact of rotational symmetry in a higher dimension.

The partial decoupling of angular momentum operators in four dimensions hints at a decomposition of the $SO\left(4\right)$ group generated by the operators L and M[48,49]. Rank 2 $SO\left(4\right)$ groups by rule can be broken down into two completely decoupled rank 1 $SU\left(2\right)$ groups: with three non-commuting generators each, yielding $3+3=6$ generators. We make then the following choice of basis and with Since both $\mathbf{I}$ and $\mathbf{K}$ constitute an $SU\left(2\right)$ known algebra, it is isomorphic to the group generated by angular momentum in three dimensions. We use our findings in Section 4.2 to get their possible eigenvalues. We can write states which are simultaneous eigenstates of H, $I_z$ and $K_z$: and one has for $i,k=0,1/2,1,\cdots$. The corresponding Casimir operators, ${\mathbf{I}}^2$ and ${\mathbf{K}}^2$, are respectively and It is important to remember that Casimir operators are distinct from generators. They are actually functions of generators and although they commute with all the generators of a group, are not considered to be in the set of commuting generators itself. A linear combination of these operators is diagonal in the basis of their common eigenvectors. The Casimir operators ${\mathbf{I}}^2$ and ${\mathbf{K}}^2$ commute with all the generators $(I_x,I_y,I_z,H,K_x,K_y,K_z)$. One has and From $\mathbf{I}$ and $\mathbf{K}$ we can construct operators and which commute with all $so\left(4\right)$ generators and H. $C_1$ yields a constraint on the eigenvalues through The last line used the relation from Equations (34) and (35). Actually, $C_1=0$ means that i and k are the same. This brings additional information (not given by $SU\left(2\right)$). We find the eigenvalues for operator $C_2$ with ease using Equations (53) and (54), yielding and Combining Equations (34), (35), (40), (51) and (52), we derive a constant expression for $C_2$: Equating Equations (58) and (59) we obtain our permissible energy levels for a hydrogenic atom: for $k=0,1/2,1,\cdots$. Recalling that $\alpha =Z{e}^2$ and setting $2k+1=n$, for $n=1,2,3,\cdots$ we recover the well-known expression for $n=1,2,3,\cdots$. The elegance of this method stems from the fact that, in addition to the hydrogen energy spectrum, the relevant constraints on various quantum numbers have emerged naturally as well. For instance, consider $\mathbf{L}=\mathbf{I}+\mathbf{K}$ from the initial definition of $\mathbf{I}$ and $\mathbf{K}$. Since $\mathbf{I}=\mathbf{K}$ , we see that by the triangular rule of vector addition, L can be at most $k+k=n-1$, and its values extend to $|k-k|=0$ in integer steps, as expected. In addition, following the algebra of $SU\left(2\right)$ groups, $I_z$ and $K_z$, analogous to $L_z$ would each have $2k+1=n$ independent eigenvalues and eigenstates. The total degeneracy of an energy level, $n\times n=n^2$, is therefore also recovered [50].

We have seen that the Runge-Lenz vector is an invariant of the Kepler motion. According the the reciprocal of the Noether theorem [1], a continuous symmetry can be associated to it. Such a symmetry was discovered by Fock and Bargmann in 1935: the Coulomb potential has an internal symmetry (or dynamical symmetry), isomorphic (for the bound states) to the one of the continuous group $SO\left(5\right)$ of the rotations in the four-dimensional space. Such a group has six parameters and there are indeed six constants of motion: the three components of the angular momentum and the three components of the Runge-Lenz vector. Actually, only five of these constants are independent, since the vectors $\mathbf{L}$ and $\mathbf{A}$ are orthogonal. The energy is not a constant of motion independent of the previous ones, since it can be calculated as a function of the invariants of $A^2$ and $L^2$, which are Casimir invariants of the $SO\left(4\right)$ group.

The harmonic potential also has a dynamical symmetry, isomorphic to the one of the continuous group $SU\left(3\right)$ of $3\times 3$ matrices with determinant equal to one. This group has eight parameters. The eight corresponding constants of motion are the three components of the angular momentum and the five components of the symmetric Laplace tensor. Five of them are independent since $\mathbf{L}$ is an eigenvector of the Laplace tensor. Let us consider a mass $\mu$ in an isotropic harmonic potential The equation of motion is It is possible to define an invariant of the motion [52] as Let us set ${\omega }^2=k/\mu$. The Laplace tensor reads [53,54,55]: where x and y are the coordinates of position operator $\mathbf{r}$ and $v_x$ and $v_y$ the coordinates of velocity $\mathbf{v}$. One has as well as yielding and and finally

We have and Eliminating v using the energy conservation gives

If both eigenvalues are positive, the orbit is an ellipse with the semi-minor and semi-major axis given by It is worth mentioning that the transformations of the group $SO\left(4\right)$ or $SU\left(3\right)$ enable one to change from an orbit of energy E to degenerated orbits of the same energy [56,57].

Analogies between the Kepler and harmonic motions were already known by Newton [58]. A mathematical transform from the one to the other has been introduced in 1965 by Kustaanheimo and Stiefel [59].

The Kustaanheimo-Stiefel transformation is a particular case of the surjective application defined by subject to the constraint $\omega =0$ with where the parameters p and t can take the values $\pm 1$. The transformation introduced by Kustaanheimo and Stiefel [60] corresponds to $p=t=-1$, while the cases $p=-t=-1$ and $p=-t=1$ (or $p=t=1$) correspond to two other (inequivalent) transformations.

The applications of the transformations (77) and (78) range from number theory to physics (classical and quantum mechanics, gauge theories, etc.). The Hurwitz matrix corresponding to (77) and (78) is connected to the problem of factoring the sum of four squared numbers. In classical mechanics, the Kustaanheimo–Stiefel transformation is used for the regularization of the Kepler problem [60]. In quantum mechanics, the latter transformation enables one to transform the Schrödinger equation for the three-dimensional hydrogen atom (in an electromagnetic field) into a Schrödinger equation for a four-dimensional isotropic harmonic oscillator (with quartic and sextic anharmonic terms) subject to a constraint [61,62,63,64,65]. Conversely, the mappings (77) and (78) may be used in some dimensional reduction process for converting a dynamical system in ${\mathbb{R}}^4$ or ${\mathbb{R}}^2\times {\mathbb{R}}^2$ into a dynamical system in ${\mathbb{R}}^3$.

It is worth mentioning that Yoshida proposed a derivation of the Kustaanheimo–Stiefel in parabolic coordinates. Through the latter derivation, it becomes clearer where and how the additional dimension is introduced and what the bilinear relation means [66].

A system with N degrees of freedom has at maximum $2N-1$ independent constants of motion. If there are N constants, the system is integrable. If there are $N+n$ ones, it is superintegrable and n times degenerated. If $n=N-1$, the system is totally degenerated, it is the case of the Keplerian and harmonic motions [67].

Prince and Eliezer have hown that the Laplace invariants are linked to the fact that the Kepler and harmonic motions are invariant by similarity, which explains respectively the third Kepler law and the isochronism of harmonic oscillations [68].

In the case of an external magnetic field, much less is known. Liberman does not provide any invariant in that case, but obtains series expansions, revealing some underlying symmetry [78]. The equation of motion of a charge q in the uniform field $\mathbf{B}$ is

Integrating with respect to time, we get where $\mathbf{k}=\mu {\mathbf{v}}_0-\mathbf{B}\wedge {\mathbf{r}}_0$ is the time-invariant Landau vector (${\mathbf{r}}_0=\mathbf{r}(t=0)$ and ${\mathbf{v}}_0=\mathbf{v}(t=0)$). Combining Equations (117) and (118), one gets which means that the motion of q is harmonic: the center C of the elliptic orbit is given by $\mathbf{OC}=(\mathbf{k}\wedge \mathbf{B})/\left(q{B}^2\right)$ and the orbit is in fact a circle since, according to Equation (117), the tangential acceleration is zero and the normal acceleration is constant. However [79], it is not necessary to know how to solve Equation (119). Multiplying Equation (118) by $\mathbf{B}$, we get where $\mathbf{v}$ is known from Equation (117). The nature of the motion is then found easily. Let us choose a new origin C with $\mathbf{OC}=\mathbf{a}$ such that the new Landau vector $\mathbf{K}$ is zero. We have, setting $\mathbf{r}=\mathbf{a}+\mathbf{R}$ and $\mathbf{v}=\mathbf{V}=d\mathbf{R}/dt$: and combining Equations (118) and (120): and where $\mathbf{L}=\mu \mathbf{r}\wedge \mathbf{v}+q{r}^2\mathbf{B}/2$ is the angular momentum $\mathbf{r}\wedge \mathbf{p}$, $\mathbf{p}=\mu \mathbf{v}+\mathbf{A}$ and $\mathbf{A}=(\mathbf{B}\wedge \mathbf{r})/2$ is the vector potential in the symmetric gauge. The invariance of $\mathbf{L}$ is a consequence of the rotational symmetry around $\mathbf{B}$ [80,81] and it is easy to show that

In the case of the motion of a charge in the field of a magnetic dipole, the motion may occur in the plane of symmetry of the dipole. If the motion is bound, the orbit fills a circular crown which axis carry the dipole: no privileged direction exists in the plane and no dynamical invariant can be found.

In the case of the motion of a charge in a magnetic monopole, there is a privileged direction in space different from the one of the cone axis: the direction of the periapsis. Such a direction is the one of a dynamical invariant which could be determined by Kerner [82,83,84]. The corresponding dynamical symmetry is the one of the $SO(4,2)$ group of matrices acting in the four-dimensional space and conserving the pseudo-norm $x^2+y^2-z^2-t^2$.

The symmetry of the quadratic Zeeman effect was investigated by different authors [76,85,86,87]. In Ref. [88] on the motion of a charge in a magnetic field, Sivardière considers the Zeeman effect for the harmonic oscillator, and shows that the motion is epicyclic.

In a very interesting work, Yuzbashyan et al., motivated by the renewed interest in non-trivial degeneracies of a simple spin Hamiltonian, proposed a general approach to extract all non-trivial symmetries from the energy spectrum [89].

The effect of a magnetic field on spectral line shapes, with or without an additional electric field, is an important topic of atomic physics (see for example the non-exhaustive list of references: [90,91,92,93,94,95,96,97,98,99]), for astrophysical or magnetic-fusion applications, both from the theoretical side and for the interpretation of Z-pinch experiments for instance.

A consistent treatment of the splitting of hydrogen atom in crossed weak electric and magnetic fields was given by Demkov et al. [100]. In that work, only the first-order perturbation theory was considered. The resulting formula for the energy corrections valid for arbitrary mutual orientation of the fields $\mathbf{F}$ and $\mathbf{B}$ was [101]: where $n^{\prime }$ and $n^{\prime \prime }$ belong to the set $\left\{-j,-j+1,\cdots ,j\right\}$, with $j=(n-1)/2$ and $e_1$ and $e_2$ are absolute values (in atomic units) of the operator: where c is the speed of light in vacuum, and The zero-order wave-functions are eigenvectors of $I_{1\alpha }$ and $I_{2\alpha }$: and where where the operators ${\mathbf{I}}_{\mathbf{1}},{\mathbf{I}}_{\mathbf{2}}$ are connected to $\mathbf{L}$ and $\mathbf{M}$ via and exactly as Equations (45) and (46). In general, first-order formula (124) completely removes the degeneracy. Second-order corrections lead only to small energy shifts. The latter was calculated by Solov’ev [102] using the dynamical symmetry group $O(4,2)$ of atomic hydrogen. The problem gets much more involved if the two fields are perpendicular ($\mathbf{F}\perp \mathbf{B}$, when $e_1=e_2=e$. The first-order correction $E^{\left(1\right)}$ depends then only on the set of quantum numbers $n^{\prime }+n=q$, implying a residual $n-$fold $(n-q)$ degeneracy. The second-order perturbation theory removes the degeneracy: where $\gamma =3ncF/H$ and $ϵ$ is the eigenvalue of the operator In the above formula, $I_{1\beta }$ is the projection of operator ${\mathbf{I}}_{\mathbf{i}}$ on a direction belonging to the plane (${\mathbf{e}}_1$, ${\mathbf{e}}_2$) and orthogonal to ${\mathbf{e}}_i$. The parameter b is given by the formula

In fact, b grows monotonously when the ration $F/H$ is increased. The extreme values $b=-3$ and $b=\infty$ correspond to purely magnetic and electric perturbations (quadratic Zeeman and Stark effects). In these limits, the quantum number q turns into the magnetic quantum number m (in the case where $F=0$) or into $n_1-n_2$, where $n_1$ and $n_2$ are parabolic quantum numbers (in the case where $H=0$).

The existence of isolated magnetic monopoles was first suggested by Dirac in 1931 [20] and could explain the quantization of the charge. Since there is still no known experimental or observational evidence that magnetic monopoles exist, they remain an active field of research. The magnetic field of the particle of charge q in $\mathbf{r}=r{\mathbf{u}}_{\mathbf{r}}$ is spherical and reads [84] and the equation of motion is where $\delta =-qg$, $\mathbf{L}=\mu \mathbf{r}\wedge \mathbf{v}$ is not an invariant, but the Poincaré vector $\mathbf{J}=\mathbf{L}+\alpha {\mathbf{u}}_r$ is one [103].