Two Masses in Relative Hooke’s Potential and Elliptic Integrals

We consider the relative motion of the system of two masses connected by a spring. We analyze it in a range of the Hooke’s law and show that the equations of the relative motion of the system are nonlinear once the equilibrium length of the spring is nonzero. Although the way of deriving the equations of motion is standard in classical mechanics solving them is a complicated and interesting problem of mathematical physics. The analysis leads naturally to elliptic integrals. We obtain complete formulas in an interesting, from both mathematical and physical point of view, way. Our analysis might be useful in some problems of molecular dynamics of diatomic molecules.


Introduction
Although the Hook's potential was discovered before the Newton's, the movement in the Hook's potential has not been thoroughly analyzed. The reason is that although it is easy to formulate the (nonlinear) system of equations describing such a motion, the solution of it is by no means easy. Essentially similar situation one encounters once analyzing a movement of a particle in the gravitational field of a black hole 1 3 Journal of Nonlinear Mathematical Physics (2022) 29:504-522 described by the Kerr'a metric. In his seminal 1973 paper [21] Teukolsky provided appropriate equations for such a movement, but at that time the solutions were obtained only by means of numerical methods and only for a movement of a particle whose velocity vector is almost tangent to the equator plane. Much later in [8] an analytical solutions for bounded orbits were given. In this paper we give a complete analytical solution of the system describing the movement of two masses in the Hooke's potential therefore filling the aforementioned gap in the literature. The difficulty of this problem, like in the above astrophysical problem is not formulation of the system of equations describing the physical phenomenon but its solution.
Both Hooke's and Newton's potentials are radial. Movement in the Newton's potential serves as a standard example of radial potential problem. In this paper we distinguish Hooke's potential U(r) = 1 2 (r − l) 2 from the potential of isotropic harmonic oscillator U(r) = 1 2 r 2 . So, the equilibrium length l of the spring is one of the parameters of the Hooke's model. The value l = 0 yields of course the potential of isotropic harmonic oscillator. We think that the problem of movement in the Hooke's potential for any l ≥ 0 was not adequately treated in the literature. In this paper we give precise mathematical formulas which describe such a motion.
In theoretical mechanics vibrating systems are of particular importance (cf. [2,9,15,17]). There is a big variety of such systems, which constitute models for various sometimes very complicated processes. Even simple one-dimensional vibrations of two masses on a spring have numerous generalizations. For example one can consider one-dimensional vibrations assuming that a spring has a finite non-zero mass (cf. [20]).
In this paper we consider two masses m 1 , m 2 connected by a spring of the equilibrium length l. The system is rotating with some angular momentum and the mass center is moving with the velocity . We assume that the spring is stiff, so that it cannot bend, and its only possible deformation is a change of length. We also assume that its mass can be neglected. We analyze this system in the range of validity of the Hooke's law. We think that this is a very interesting system also from the theoretical point of view. By a standard procedure applicable for a radial potential one can separate the total motion of the system into the motion of the mass center and the relative motion. This is recalled in Sect. 2.
The relative motion takes place in the plane perpendicular to the vector of relative angular momentum, whereas the velocity vector of the mass center motion is arbitrary (cf. Fig. 1). Naturally, in our analysis we use the basic conservations laws of energy and angular momentum. It is well known that systems with any radial potential are integrable, but the quadratures are not expressible by elementary integrals except for few forms of potentials. We show here that even in the Hooke's law range the quadratures are not elementary. Our approach allows one to consider the relative motion of this system as a collective effect of rotation and vibration. Qualitatively, the movement has completely different form than that analyzed under a "small oscillations" assumption in a neigbourhood of the minimal potential. This can be justified by the fact that the Hessian has singularities in the points which render its minima. This corresponds to the radial symmetry in the plane of relative motion (p.r.m.). Minimal value of relative potential is zero and it is attained on the circle r = l in p.r.m. On the other hand the "small oscillations" assumption imply that the Hessian has isolated critical points at which it is nonsingular.
Since the spring in the equilibrium state has a nonzero length the equations describing rotation and vibrations are nonlinear, even in the range of validity of the Hooke's law. Therefore the dynamics of both motions is inseparable in the sense of linear superposition although the geometric separation of these motions is clear. It is worth noting, that the above mentioned collective effect leads to much more elaborated formulas than those for each of the motions (rotation and vibration) alone. We provide exact formulas for this in terms of elliptic integrals. The formulas were programmed and the results graphically depicted by use of [23]. We suggest a solution of differential equations of the movement in a parametric form. Such an approach was both used in the solution of the classical problems as the Kepler's problem, and contemporarily, for describing motion of a particle in the Kerr's metric (cf. [8]). However, this approach is fruitful, as we show in our paper, for analysis of rotating vibrational systems. In this paper we take into account mutual influences of rotational and vibrational components in a relative motion for a relative Hooke's potential. This type of analysis is necessary for the description of the dynamics of the system. From theoretical point of view, we think, this is good approach to analysis of dynamics of a diatomic molecule. Our solution of the model is exact. From mathematical point of view the system is interesting as the relative motion is described by means of elliptic integrals of the first and third kinds. The trajectory of motion can be viewed as a generalized Lissajous curve (cf. [12]). In the analysis one can use the tables of elliptic integrals [5] (cf. also [1,11]). We derive appropriate formulas by means of algebraic and trigonometric transformations. The advantage of our approach is that the formulas we obtained guarantee no computational problems once one changes branches of square roots. This is often the case when one uses formulas from tables directly. The content of the paper is the following. In Sects. 2 and 3 we recall basic well known facts concerning the movement of two masses in the radial potential. In Sect. 4 we derive the system of differential equations describing the movement of our system. We define there certain parameter Δ whose sign governs the movement both qualitatively and quantitatively. The heart of the paper are Sects. 5-7. Each of these sections are technical. They start with the choice of a parameter used to obtain the parametric solution for the considered case and then we show how the solution is derived. In Sect. 8 we include additional comments on our results and briefly discuss some nonlinear perturbations of the Hooke's law.

Separation of Motion into Mass Center and Relative Motions
Let us recall the standard procedure of separation of motion into the motions indicated in the title of this section. The Lagrange function is given by the following formula: where is a radial potential.
The mass center motion satisfies the condition of inertia law. The relative motion is equivalent to the motion in a central radial field. Indeed, it is very well known that the separation of the motion into the mass center and relative parts leads to the following decomposition: where L mc ( L r resp.) is the mass center part of the Lagrange function (relative part of Lagrange function resp.). Similarly, mc ( resp.) is the mass center position (relative position resp.). Quantity m is the reduced mass of m 1 and m 2 .

General Case of Motion in Radial Potential
In the radial potential field the equations of motion of a material point are given by the following nonlinear system of ordinary differential equations [15]: obtained from the conservation laws of energy and angular momentum Notice that in this case the effective radial potential U ef f equals J 2 mr 2 + U(r). Hence where Let the condition be such that the radial motion is limited to r ∈ ⟨r 1 , r 2 ⟩ , where r 1 , r 2 are zeroes of Q(r):

This yields
Notice that once the range of the motion is known r ∈< r 1 , r 2 > the energy and the angular momentum are determined by the formula (13).
In the sequel we use elliptic integrals of the first and third kinds. We adopt the notation from [23]. Notice that it differs from the convention used in [5]. Thus is the (incomplete) Legendre elliptic integral of the first kind and is the (incomplete) Legendre elliptic integral of the third kind.

Spring in a Harmonic Range
In this section we obtain the system of nonlinear differential equations for the motion in the Hooke's potential. In the following sections we solve them parametrically.
For the harmonic case one gets Equivalent form of the formula for energy is the following Now, let us determine the remaining roots r 3 , r 4 of the quartic (17). The Vieta formulas and (18) lead to following formulas: The roots r 3 , r 4 are given by the following formula: where The differential equations of motion can be written as the following equalities of differential forms: where and a = r Remark 1 Notice that our solutions for various signs of Δ are obtained by careful choices of the parameters. Although in each case the parameter is different we denote them by the common letter .  Moreover, the formula for r is given by (29). In formulas (36) and (37) the elliptic integral of the first kind is given by (14) and the elliptic integral of the third kind is given by (15).
In Figs. 2, 3 and 4 we depicted functions t( ), r( ) and ( ) computed, for a particular data, from formulas (36), (29) and (37). Since for other cases eg. Δ = 0 and Δ < 0 the character of these functions (not the exact values) is similar we do not show them there. In Fig. 5 (Fig. 6 resp.) we displayed trajectories in the x − y plane for the same data and the range of −10 < < 10 ( −12 < < 12 resp.). This allows us to see how the trajectory comes into existence.
. Notice that from formula (24) we see that Δ = 0 if and only if r 1 + r 2 − 2l = 0 or (r 2 −r 1 ) 2 r 1 +r 2 − 2l = 0. The first equality leads to radial harmonic oscillations without rotation. We omit this case as a trivial one. For . . (43) (45) where K l.c. (x) = x − arc tan(tan x) is a locally constant function with such a jump that it is possible to make I p a continuous function (by assigning appropriate values in the points which are not in the domain). So and the parametrized equations of the motion one gets from (25) : and r( ) is given by (40). An example of a trajectory for Δ = 0 is depicted in Fig. 7. Notice that in Fig. 7 we see white spots which correspond to points with removable discontinuities cf. formulas (45) and (47).

The Case 1 < 0
For the case Δ < 0 the roots c, d are complex and conjugate. Hence, we consider the following form (53)

Parametrized Solution of Motion Equations for 1 < 0
From the formulas of the previous subsection in particular from (52) and (55) (cf. also (26)) finally one gets the parametric formulas for the relative motion. The parametrization of time is the following: The parametrization of is given by the following formula: The function r( ) is given by (50).
The function in the last summands of (57, 58) f 1 (H, , k 2 ) is an elementary integral of the following type: By means of the Abel substitution we obtain Hence (56) In Fig. 8 we showed a trajectory for Δ < 0 and an exemplary data.

Spirograph-Like Orbits
Notice that in the x − y coordinates in the p.r.m we have obtained spirograph-like orbits. From an old theorem of Bertrand [3] it follows that for the radial potential U = kr and the bounded motion r min ≤ r ≤ r max the closed orbits in the phase space are obtained only for = −1 and = 2 . In the first case we obtain the force corresponding to gravitational or Coulomb law in the second case the force obeys the Hooke's law, see [10,15] and especially [19]. In [19] the authors also consider a circular motion of radius r 0 in the radial potential U = kr along with the radial motion in the Hooke's potential of a harmonic oscillator i.e.
where k ef f = ( + 3) J 2 mr 4 0 and J is a constant angular momentum.
This in fact corresponds to the situation of linear perturbation of force. Then the authors separate both motions (vibration and rotation) and compute the ratio of the periods of these motions. If this ratio is a rational number then the orbit is closed.
Geometrically this situation corresponds to the well known situation of winding a half-line with constant speed on a two-dimensional torus which rotates with the constant angular speed about it axes of symmetry (cf. [2]) Remark 2 Notice that if for the aforementioned separation of motions the orbit is closed in the phase space then it is obviously closed in the x − y plane.
System we consider is much more general, we do not make any limitations (except that the spring is massless) on its parameters. Especially, we do not assume that the motions are separable. In fact they are not. In Fig. 6. we depicted the situation where the angular momentum is relatively small but the orbits are non-closed. One can see the "shift" of the orbit near the x axis.The non-closure of the orbits is a typical situation for our potential. As it is mentioned in [19] closed orbits are possible only for a discrete set of parameters. In Figs. 7 and 8 we depicted non-closed orbits for Δ = 0 and Δ < 0. The character of the bounded orbits for the potential U = kr is the same for > −3 . (For = −3 one may obtain a Cotes' spiral as an orbit) [19]. However, since the parameterizations for Δ > 0 , Δ = 0 and Δ < 0 differ essentially, the exact coordinates and dynamics for such orbits are different. Our solution gives the exact values of the coordinates.

Some Related Physical Systems
It is well known that a simple pendulum only for small angles of a swing can be modeled by the equation of a harmonic oscillator. This is justified by the well known limit: lim →0 sin = 1. However, once is not close to zero (one cannot assume any longer that sin ≈ for ≤ ) and the solution of the equation of motion is given by the following formula (cf. [16]): 1 − sin 2 1 2 sin 2 = sn −1 sin , sin 1 2 In the above formula is an angle between the vertical and the pendulum, is the amplitude of the angle of swing, l denotes the length of a pendulum and g is the gravitational acceleration. Another one dimensional system, where the Hooke's potential is disturbed by the quartic term, is the potential of the Duffing oscillator. Then the tension in the spring is given by the following equation: where x is the extension of the spring. If c > 0 then the spring is called stiff whereas for c < 0 it is called soft. The equation of motion, after appropriate rescaling of time, can be put into the following form [16]: The solution for the motion of the stiff spring with the initial conditions: x(0) = a ≠ 0 and ̇x(0) = 0 is given by the following formula involving the elliptic cosine function [16]: For the soft spring the solution can be expressed by the elliptic sine function as follows: The above described two (pendulum and Duffing oscillator) 1-dimensional systems with the perturbations of the Hooke's potential require elliptic functions, directly related to elliptic integrals, for expressing the motion of a mass. For a generalization of the simple pendulum see e.g. [13] and for generalizations of the Duffing oscillator see [6,7,14].
Notice that in our system we use the Hooke's potential in the form U(r) = 1 2 (r − l) 2 , l > 0 rather than the potential of the harmonic oscillator (the case l = 0 ). When the system only vibrates this distinction corresponds to changing the coordinate on the line and therefore is irrelevant. However, once the angular momentum J is nonzero we obtain the 2-dimensional nonlinear system given by the coupled equations (6) and (7). We have not seen analysis of such a system in the literature. We obtained parametric solutions for each of the following cases: Δ > 0 , Δ = 0 and Δ < 0, where Δ is given by the formula (24). One can try to solve the system of differential equations using one of the finite difference schemes. In fact checking our formulas we used the Runge-Kutta method of order 4 for the coupled initial value problem (see [4]). The differences of the results were negligible. It should be noted that as usual the difference   method will be working well for not too big intervals [t min , t max ] of time. Our method has no such limitations. We obtained exact formulas for any value of the parameter .

Conclusion
In the paper we have analyzed a system of two masses connected by a massless spring. We took into account both rotation and vibration of the system. In our opinion the research presented here might be useful in some problems of molecular dynamics of diatomic molecules (cf. [18,22]). In such models the assumption that the spring is massless is natural from the point of view of classical theory. For polyatomic molecules a typical way of taking into account mutual influences of rotation and vibration is a semi-rigid molecule model. These influences (in a suitable reference frame ) are given by a Coriolis coupling and varying with time moment of inertia (cf. [22]).This model can be solved numerically producing satisfactory results for small amplitudes of oscillations. Analyzing the Coriolis coupling for diatomic molecules we readily see that it vanishes. In this case the rotational-vibrational coupling is given by change of the moment of inertia.
Our solution is exact for the model under consideration. Notice that although we used, in our calculation, the same letter for an independent parameter of motion, there is a discontinuous change of when the sign of the discriminant Δ is being changed. However, this does not effect the character of trajectories. This is a significant difference in comparison to Kepler's problem, where the change of sign of the discriminant results in a transition form a bounded trajectory to an unbounded one (or vice versa).