Article The Free Euler Rigid Body Revisited Lidia Jiménez-Lara 1 and Jaume Llibre 2, * 1 Departamento de Física, Universidad Autónoma Metropolitana–Iztapalapa, P.O. Box 55-534, Mexico D.F. 09340, Mexico; lidia@xanum.uam.mx 2 Departament de Matemàtiques, Universitat Autònoma de Barcelona, 08193 Bellaterra, Spain * Correspondence: jaume.llibre@uab.cat Abstract: We review from a different perspective the approach and solution to the torque-free Euler equations, also called the free asymmetric top equations. We aim to simplify and broaden the study of the asymmetric free rigid body. This is an old but important integrable problem that has two first integrals: the energy and the angular momentum. We reduce this problem by eliminating the time as the independent variable in the three autonomous Euler equations written in cylindrical dimensionless variables, which allows a geometric study of the solution as a function of the cylindrical angle variable ψ, by means of continuous deformations dependent on the two independent parameters κ and e0 . The parameter space is divided into six disjoint regions, whose boundaries are the separatices and degenerated cases. The solutions are given in terms of trigonometric functions of the independent cylindric angle ψ. Keywords: free rigid body; Euler top; reduction; energy hyperboloid; global solution 1. Introduction The study of the rigid body dates back to the XVIII century, with the pioneering works of Euler (1707–1783), who dedicated a large part of his life to study this problem, trying to explain and predict the motion of ships and their building (a historical review of Euler’s works on the rigid body are in the article by Marquina et al. [1]). Around 1736, Euler published in two large volumes his treatise Mechanica sive motus scientia analytice exposita (Mechanics or the science of motion, expounded analytically), where he proposes that Citation: Jiménez-Lara, L.; Llibre, J. the motion of a rigid body can be studied as two types of combined movements: one of The Free Euler Rigid Body Revisited. translation around its center of gravity, and another of rotation of the body’s orientation Dynamics 2023, 3, 603–621. https:// around an axis that passes through its center of gravity. The most simple rotation problem doi.org/10.3390/dynamics3030032 is when there are no torques, the so called Euler free rigid body or simply the Euler top. In 1755, Segner showed that every rigid body has three main axes where the inertial tensor Academic Editor: Christos Volos is diagonal. Then, in the series of works that Euler wrote between 1758 and 1765, which Received: 4 August 2023 culminate in the voluminous compendium Theoria motus corporum rigorumum seu solidorum, Revised: 8 September 2023 Euler used the principal axes system, and found the so-called Euler equations and the Accepted: 12 September 2023 Euler angles. He obtained the kinematic relations for the motion of a heavy rigid body Published: 18 September 2023 and addressed the integrable problem of the free rigid body, reducing it to quadratures. The angular momentum in the inertial frame is a constant of motion, but in the body-fix frame it satisfies the Euler equations. A geometric construction of the solution to the Euler equations for the free rigid body Copyright: © 2023 by the authors. was given by Luis Poinsot [2] in 1851 using the polhod and herpolhod: the inertia ellipsoid Licensee MDPI, Basel, Switzerland. This article is an open access article with a fixed center rolls without slipping on a fixed plane perpendicular to the constant distributed under the terms and angular momentum; the position vector of the contact point is like an instantaneous axis of conditions of the Creative Commons rotation, and draws polhodes on the ellipsoid of inertia, or herpolhodes on the fixed plane. Attribution (CC BY) license (https:// The general solution of the motion of a rigid body around a fixed point under no creativecommons.org/licenses/by/ forces was probably first given by Rueb in 1834 [3] in his PhD Thesis, and completed by 4.0/). Jacobi in 1849 [4], who obtained the Euler angles in terms of elliptic functions and theta Dynamics 2023, 3, 603–621. https://doi.org/10.3390/dynamics3030032 https://www.mdpi.com/journal/dynamics Dynamics 2023, 3 604 functions, and their expression in terms of Fourier series. Elliptic integrals had already appeared in the Elastica of Jakob Bernoulli, while Maclaurin, Fragano, Legendre and others used them in connection with the problem of rectifying an arc of an ellipse. Proceeding in analogy with the circular functions, Abel, Jacobi, Gauss and Weierstrass [5] introduced the idea of inverting the elliptic integrals defining the elliptic functions. Carl Gauss noticed that the rigid body position may be uniquely determined by quaternions with unit norm, now known as the Rodrigues–Hamilton parameters [6]. Integrable cases of the rigid body were the first to be studied. Besides the free Euler top, just two integrable cases for the Euler equations with an additional global real integral are known: the Lagrange case (1788) [7], a symmetrical top with the center of mass on the symmetry axis moving around a fixed point in a constant gravitational field, and the Kovalevskaya case (1888) [8], a symmetrical top with I1 = I2 = 2I3 and the center of mass lying in the plane perpendicular to the symmetry axis. The Goryachev case (1900) [9], a symmetrical top with I1 = I2 = 4I3 and the center of mass lying in the equatorial plane, has only particular integrals. It was integrated by Chaplygin in terms of hyperelliptic integrals [10]. The Goryachev–Chaplygin case has a rational relationship with the three- body periodic Toda lattice [11]. The additional integrals in the Euler and the Lagrange cases are related with natural physical quantities: the square of the angular momentum, and the projection of the angular momentum on the dynamical symmetry axis, respectively. The additional integral in the Kovalevskaya case is highly non-trivial: it was found a hundred years later than Lagrange’s case. The basic methods for finding first integrals from the symmetries and for studying the integrability are the separation of variables and the Noether’s theorem. A new qualitative frame to study integrability in finite dimensional Hamiltonian systems was shown by Arnold [12,13]. In his course on Celestial Mechanics, Andoyer (1926) [14] used the two integrals of motion for the Euler and Lagrange integrable cases (the angular momentum and its projection into the z body-axis) as convenient variables to study the rotational motions of the planets. Similar variables were introduced by Deprit [15] in 1967, reducing the free rotor system to a conservative Hamiltonian system with one degree of freedom. The phase space has a pendulum-type dynamics. These variables are widely used in the perturbation theory of the rigid body. Lagrangian and Hamiltonian structure of the rigid bodies can be studied with Euler angles, for instance, but many questions are easier to solve using the three Euler equations of motion for the angular velocity or the angular momentum. The dynamics of the rigid body has played a major role in the development of modern geometrical methods and theories of Lie groups and algebras [16,17]. The rigid body equations can be written in Hamiltonian form as a Lie-Poisson system associated to a Lie algebra structure in R3 . Differ- ent algebras can be chosen: SO(3), SE(2), SE(2, 1). The deformation through these algebras can be achieved defining SL(2, R)-linear combinations of the constants of motion (Casimir constants) [18,19], the energy and square of the body angular momentum. Under these linear transformations, the Euler equations of motion and the trajectories in R3 remain unchanged. Choosing the SL(2, R)-linear transformation properly [18], the dynamics of a simple pendulum is obtained. Separability in the Liouville sense became an important issue with the advent of quantum mechanics. Reiche [20] opens a research path in this direction, showing that the classical kinetic energy of the free asymmetric Euler top is separable whenever the canonical momentum corresponding to the gyration angle around the direction of the constant angular momentum is zero; the intersection of the ellipsoid of energy E and the sphere of angular momentum L2 can be parametrized in sphero-conical coordinates in terms of Jacobi elliptic functions. In these coordinates, the asymmetric free rigid body can be quantized, and the Schrödinger equation separated [21–26]. For more recent works on the rigid body see [27–32]. In this work, we study the trajectories of the solutions by eliminating the time as the independent variable in the autonomous Euler equations. It is convenient to use the di- Dynamics 2023, 3 605 mensionless variables and independent parameters of asymmetry and energy-momentum, as in the references [5,19,33], as these allow for a geometric study by means of continuous variations on the two independent parameters. We intended this work to be self-contained. In Section 2, we write the Euler equations in dimensionless variables and dimensionless parameters, the asymmetry parameter κ, and the energy-momentum parameter e0 . In these variables, the “energy surface” becomes generically a hyperboloid. The geometry of the hyperboloids depends on both parameters κ and e0 . The asymmetry parameter κ fully specifies the values of the three dimensionless inertia parameters e1 , e2 , e3 , which divide the space of parameters into six regions; for a given κ, the geometry of the hyperboloid depends on the relative values of e0 and the inertia parameters, as expected. Special cases take place when the value of e0 is on the boundary of these regions, i.e., e0 = ei for some i = 1, 2, 3. In Section 3, we use cylindrical coordinates and eliminate the time as an independent variable in the autonomous Euler equations to obtain the equation of the trajectory in a two-dimensional configuration space, with non-autonomous differential equations, which are solved explicitly in terms of trigonometric functions of the cylindric angle ψ. Of course, the solutions depend on where in the six regions, or their boundaries, the values of e0 and κ parameters lie. We give the time-independent trajectory solution as a function of the two continuous parameters. Continuously changing the values of κ and e0 , we obtain all the solutions for the trajectory and the geometry of the hyperboloid of energy for all the possible cases of symmetries and asymmetries of the rigid body. In Section 4, the solutions are studied as limit cases of Section 3, at the boundary of the six open regions. They are the separatrices of different kinds of solutions, which are changing their symmetry axis. In Section 5, the conclusions are summarized. 2. The Euler Equations To study the rotational motion of a torque free rigid body, we assume that the body’s center of mass is fixed at the origin of an inertial frame. The Newton equations are dJ T= , (1) dt where T is the torque, and J the angular momentum of the rigid body in the inertial frame. If T = 0, J becomes a constant of motion. These vectors in the body-fix reference frame are denoted by N and L, respectively, and they are related by an orthogonal rotation matrix R, such that T = R N and J = R L. Writing the Equation (1) in the body frame, we have the Euler equations N = ω × L + L̇, (2) where ω is the pseudovector angular velocity, such that Ṙ a = R ω × a for all vectors a in the body frame. This equation gives the velocity due to the rotation of a point in a. If the torque vanishes, either because the motion is free of external forces, or be- cause the rigid body is in a constant field and the origin of the non-inertial frame is the body’s center of mass, then T = N = 0 and we are dealing with the rigid body torque free system with two integrals of motion: the energy E and the angular momentum J (then J = L = constant). In this case, the free Euler equations are given by L̇ = L × ω. (3) Dynamics 2023, 3 606 Since L = Iω with I the matrix of moment of inertia, the Euler equations for the torque free rigid body can be written either as L̇ = L × I−1 L. (4) or I ω̇ = I ω × ω. For tensor I to be diagonal, we choose the axes of the body frame in the direction of the principal axes system. The special case with spherical symmetry I1 = I2 = I3 has the solution L = const., and the rotation matrix R is also constant. The kinetic energy E and the square of the angular momentum L2 are given in terms of the moments of inertia and the components of L: 1 T 1 L21 L2 L2 E= ω I ω = L T I−1 L = + 2 + 3, (5) 2 2 2I1 2I2 2I3 L2 = L2 = L21 + L22 + L23 . (6) Geometrically, the equation √ for the conservation √ √ of the energy generically defines an ellip- soid with semi-axes 2EI1 , 2EI2 and 2EI3 , whereas the equation for the conservation of the modulus of the angular momentum defines a sphere of radius || L|| = L. If a solution exists, the radius of the sphere must have a value between the minimum and maximum values of the semi-axes of the ellipsoid 1 2E 1 ≥ 2 ≥ . (7) Imin L Imax When the vector L moves relative to the top’s principal axes of inertia, it lies along the curve of the intersection of the surfaces E = const. and L2 = const. Since the paths are closed, the motion of L with respect to the body system must be periodic, describing some conical surface and returning to its original position. To introduce the dimensionless variables [19,33], let u be in the unitary vector in the direction of the angular momentum L, with the equation of motion u̇ = L u × I−1 u. (8) If the tensor I−1 = diag{1/I1 , 1/I2 , 1/I3 } is written in terms of its irreducible representa- tions (its trace and a traceless symmetric matrix) [19] 3   3 1 1 1 L T I−1 L = Tr(I−1 ) L2 + ∑ − Tr(I−1 ) L2i = TL2 + A ∑ ei L2i , (9) 3 i =1 Ii 3 i =1 where A and T have dimensions of the inverse of the momentum of inertia 1 1 T= Tr(I−1 ), − T = A ei , i = 1, 2, 3, (10) 3 Ii the inertia parameters ei here defined are dimensionless. We find that 3 1 1 1 ∑ A ei = I1 + + − 3T = 0 I2 I3 i =1 and then the trace of the dimensionless matrix   e1 0 0 E =  0 e2 0 0 0 e3 Dynamics 2023, 3 607 is canceled e1 + e2 + e3 = 0 . Computing the trace of the square of the matrix AE, we have " # 2 1 1 1 1 1 1 A2 (e12 + e22 + e32 ) = + 2+ 2− − − 3 I12 I2 I3 I1 I2 I1 I3 I2 I3 "  # 1 2 1 2 1 2     1 1 1 1 = − + − + − (11) 3 I1 I2 I2 I3 I3 I1 and if we take 3 e12 + e22 + e32 = 2 then " # 4 1 1 1 1 1 1 A2 = + 2+ 2− − − . 9 I12 I2 I3 I1 I2 I1 I3 I2 I3 In order that A 6= 0 and ei be well defined, from (11) we see that the rigid body cannot be completely symmetric; in what follows we exclude this case, but the symmetric case with only two equal moments of inertia is included in our discussion. The inertia parameters turn out to be convenient because they are bounded and satisfy the equations e1 + e2 + e3 = 0, 3 (12) e12 + e22 + e32 = . 2 In the three dimensional space of inertia parameters, these equations define√ a plane going through the origin, and a sphere with its center at the origin and radius 3/2, respectively. The main difference between the dimensionless inertia parameters ei and the inverse moments of inertia Ii−1 is that the ei parameters can be positive, negative or zero, while the Ii−1 are always positive reals. However, the equations of motion and their solutions rely only on the differences of the inertia parameters ei − e j and, subsequently, on the differences 1/Ii − 1/Ij as per Equation (10). Furthermore, since the sum of the inertia parameters must be zero, two of them must possess opposite signs, while the third one can have any sign, including perhaps being zero. Notice that if one inertia parameter is smaller than another, like ei < e j , it means that their moments of inertia are in the opposite order Ii > Ij , because ei is related to 1/Ii . According to the conditions (12), only one parameter is needed to determine the three inertia parameters, named the asymmetry parameter κ such that e1 = cos(κ ) , e2 = cos(κ − 2π/3), e3 = cos(κ + 2π/3). (13) The angle κ parametrizes the circle defined by the intersection of the plane and the sphere in the Equation (12). A third equation also holds, e1 e2 e3 = 14 cos(3κ ), related with the third invariant (the determinant) of the characteristic polynomial of the inverse inertia matrix, see [19] for more details. The dimensionless parameter e0 associated with the energy and the square of the angular momentum is defined analogously to (10) 2E − T = Ae0 , (14) L2 and is called the dimensionless energy-momentum parameter or shortly the energy parameter, although it depends on both energy and the squared of the angular momentum. Dynamics 2023, 3 608 The conservation of energy (5) and angular momentum (6) in the new dimensionless variables have been transformed into e1 u21 + e2 u22 + e3 u23 = e0 (15) u21 + u22 + u23 = 1 (16)  and the Euler Equation (8) results in u̇i = LAε ijk ek − e j u j uk . Finally, the dimensionless time is defined by τ ≡ tLA in order to obtain the equations dui  = ε ijk ek − e j u j uk . dτ From now on, a dot will denote the derivative with respect to the dimensionless time τ u̇i = ε ijk (ek − e j )u j uk . (17) These are the dimensionless Euler equations: all the variables ui , “time” τ and inertia parameters ei are dimensionless, as well as the energy-momentum e0 , Equation (14). The surface (15) does not correspond to the energy ellipsoid (5), but they are related through the transformation (14). This is why in the new parameters e1 , e2 , e3 that satisfy (12) do not define an ellipsoid but a hyperboloid. Abusing the language, we call e0 the energy parameter (which also involves the square of the angular momentum), and the equation e1 u21 + e2 u22 + e3 u23 = e0 is called the energy level. For any value of κ, at least one of the inertia parameters is positive and another is negative emin < 0 < emax . The third one may be either positive, negative or zero. For instance, if κ ∈ [0, π/3], the inertia parameters satisfy e3 ≤ e2 ≤ e1 (with e3 < 0, e1 > 0, and e2 either positive, negative or zero), and the moments of inertia are in the opposite ordering I1 ≤ I2 ≤ I3 , with I1 6= I3 because we assume that the rigid body has no spherical symmetry. In particular, when the asymmetry parameter is a multiple of π/3, mod(2π)   π 2π 4π 5π κ ∈ Ks = 0, , , π, , , (18) 3 3 3 3 some ei have an extremal value and the rigid body is symmetric, with two equal inertia parameters. Likewise, when the asymmetry parameter is an odd multiple of π/6, mod(2π)   π π 5π 7π 3π 11π κ ∈ Ka = , , , , , , (19) 6 2 6 6 2 6 the rigid body is the most asymmetric one: one inertia parameter is zero, and the other two inertia parameters have the same absolute value but opposite signs. Figure 1 shows the inertia parameters as a function of κ. These special values of κ split the interval [0, 2π ) in twelve disjoint open subintervals. For the values of κ ∈ [0, 2π ), we can continuously cover all the possible asymmetries of the rigid body: symmetrical top (κ ∈ Ks ) and asymmetrical top (κ ∈ / Ks , including the most asymmetrical tops κ ∈ Ka ), with all the possible ordering of the ei ’s parameters, as it is shown in the Figure 1; of course, instead of the inverse inertia matrix I−1 , we are using the dimensionless matrix E = diag{e1 , e2 , e3 }. The transition from three real inertia moments to a single asymmetry parameter and the conversion of the ratio E/L2 ∈ R into the energy-momentum parameter e0 have allowed us to reduce all case studies to only two bounded independent parameters, namely κ and e0 . For κ taking continuous values within its domain [0, 2π ), we naturally obtain the six possible permutations of the three inertia parameters, ranging from symmetric rigid bodies (treated as subsections (c) in all the studied cases, where κ is a multiple of π/3), and notably, to the most asymmetric cases (when one inertia parameter vanishes, if κ is an odd multiple of π/6). Dynamics 2023, 3 609 1 e1 e2 e3 2ic 3iic 1iic 0.5 1iib 2ia 2ib 3iia 3iib 1iia π π π 2π 5π 7π 4π 3π 5π 11 π π 2π 6 3 2 3 6 6 3 2 3 6 -0.5 3ia 3ib 1ia 1ib 2iia 2iib 3ic 1ic 2iic -1 Figure 1. Parameter space showing Case 1 with solution (42), Case 2 with solution (55), and Case 3 with solution (66). The red line represents the value of e1 , the green line represents e2 , and the blue line represents e3 . Special cases take place when the rigid body has cylindric symmetry, κ ∈ Ks . If ei 6= e j = ek for i 6= j 6= k in the set {1, 2, 3}, from (17) we obtain ui as constant and u̇ j = α uk , u̇k = −α u j , where α = (ei − e j )ui , and the solution is given in terms of circular functions. The surface of energy (15) and angular momentum (16) define the level surfaces E = {u ∈ R3 |u T E u = e1 u21 + e2 u22 + e3 u23 = e0 }, (20) S = {u ∈ R3 |u T u = u21 + u22 + u23 = 1}. (21) The sphere of angular momentum (21) is a unitary sphere, and the level surface (20), is generically an elliptic hyperboloid instead of the original ellipsoid (5), because at least one inertia parameter is positive and another is negative. Remember that the name energy surface is used for simplicity, and e0 depends also on the square of the angular momentum, see Equation (14). The solution to the dimensionless Euler equation lives at the intersection E ∩ S. The condition for the intersection between the sphere of angular momentum and the original ellipsoid of energy (7) becomes in these parameters the intersection of the hyperboloid (20) and the unitary sphere (21), with the condition emin ≤ e0 ≤ emax . (22) When e0 = emin or e0 = emax , the hyperboloid and the unitary sphere are tangent at only two points. The geometry of the energy surface (20) depends on the values of e0 and κ. We have two cases: I. e0 6= 0. The energy surface E has the equation e1 2 e2 2 e3 2 u + u + u = 1, (23) e0 1 e0 2 e0 3 which is either a hyperboloid or a cylinder, depending on the value of the asymme- try parameter: Dynamics 2023, 3 610 (i) κ∈ / Ks ∪ Ka . This is the generic case. All the coefficients in (23) are not zero, and - If two coefficients are positive, and the other one is negative, let us say ei /e0 < 0, then E is an elliptic hyperboloid of one sheet around the axis ui . - If two coefficients are negative, and the other one is positive, let us say ei /e0 > 0, then E is an elliptic hyperboloid of two sheets around the axis ui . That two coefficients in (23) are positive or negative depends on the sign e0 and on the intermediate inertia parameter, because the other two always have opposite signs. (ii) κ ∈ Ks and two inertia parameters are equal. All the coefficients in (23) are not zero, and - If two coefficients have the same positive value, and the other one is negative, let us say ei /e0 < 0, then E is a hyperboloid of revolution of one sheet around the axis ui . - If two coefficients have the same negative value, and the other one is positive, let us say ei /e0 > 0, then E is a hyperboloid of revolution of two sheets around the axis ui . (iii) κ ∈ Ka . A coefficient is zero, let us say ei /e0 = 0, and the other two have the same absolute value but opposite signs. This is a hyperbolic cylinder with cylindrical axis ui . II. e0 = 0. The energy surface E has the equation e1 u21 + e2 u22 + e3 u23 = 0, (24) which is either a cone or two intersecting planes, depending on the value of the asymmetry parameter: (i) κ ∈/ Ks ∪ Ka . All the coefficients in (24) are not zero. If two coefficients have the same sign, and the other the opposite one, let us say ei . E is an elliptic cone around the axis ui . (ii) κ ∈ Ks . All the coefficients in (24) are not zero: two coefficients have the same value and sign, and the other has the opposite sign, let us say ei , then E is a circular cone around the axis ui . (iii) κ ∈ Ka . Then some inertia parameter vanishes, let us say ei = 0, and E are two planes intersecting along the ui axis. The solution of (17) for asymmetric rigid bodies is well known in terms of elliptic Jacobi functions, see for instance [33,34], and they live at the intersection of the unitary sphere (6) and the hyperboloid of energy (5). In the next section, we will take a different approach. We will explicitly find the trajectory equation as the solution to the Euler equations in a reduced two-dimensional space, where the solution depends on trigonometric functions of a cylindric angle ψ. The geometry of the solutions is related with the hyperboloids previously characterized in terms of the asymmetry parameter κ and the energy parameter e0 . 3. The Explicit Solution for the Trajectory Equation We write the dimensionless Euler Equation (17) in the abbreviated form u̇1 = ( e3 − e2 ) u 2 u 3 = a u 2 u 3 , u̇2 = ( e1 − e3 ) u 1 u 3 = b u 1 u 3 , (25) u̇3 = ( e2 − e1 ) u 1 u 2 = c u 1 u 2 , Dynamics 2023, 3 611 where the coefficients a, b and c are the differences between the inertia parameters, which satisfies a + b + c = 0 and can be written in different ways √ √ a = e3 − e2 = − 3 sin(κ ) = 3 cos(κ + π/2), (26) 1 √  √ b = e1 − e3 = 3 sin κ + 3 cos κ = 3 cos(κ − π/6), (27) 2 1 √  √ c = e2 − e1 = 3 sin κ − 3 cos κ = − 3 cos(κ + π/6). (28) 2 We expect these differences of inertia parameters to be present in the solution. We will introduce cylindrical coordinates for distinguish three possible cases depend- ing on the symmetry axis of the cylindrical coordinates. Then we will solve the Euler Equation (17) eliminating the time τ to find the trajectory equation parametrized with a cylindric angle ψ. In what follows, we assume r ≥ 0. In Figure 1, a summary of all the cases and subcases is shown. Case 1. Cylindric coordinates (u1 , r, ψ) around the u1 axis u1 = u1 , u2 = r cos ψ, u3 = r sin ψ. (29) Equation (25) become u2 u˙2 + u3 u˙3 ṙ = = (b + c) r u1 sin ψ cos ψ, r u2 u˙3 − u3 u˙2 (30) = u1 c cos2 ψ − b sin2 ψ ,  ψ̇ = r2 u˙1 = a r2 sin ψ cos ψ. Note that u̇1 ṙ = − u , i.e., r ṙ + u1 u̇1 = 0, r 1 and we have the first integral u21 + r2 = const., where r2 = u22 + u23 , which is the conservation of the norm of the angular momentum given in (21) u21 + r2 = 1. (31) To solve the equations of motion (30), we change to the independent variable ψ, reducing the differential system in one dimension dr (ψ) ṙ r (ψ)(b + c) sin ψ cos ψ r 0 (ψ) = = = , (32) dψ ψ̇ c cos2 ψ − b sin2 ψ du1 (ψ) u˙ ar2 (ψ) sin ψ cos ψ u10 (ψ) = = 1 = . (33) dψ ψ̇ u1 (ψ) c cos2 ψ − b sin2 ψ These are the trajectory differential equations in the independent variable ψ. Now the system is non-autonomous in a reduced two-dimensional space (r, u1 ). In order that these equations be well defined, the denominators of (32) and (33) must be different from zero for all ψ. We assume that u1 (ψ) 6= 0 and b c < 0, (34) and also r < 1 by (31). Critical points take place when r = 0 at the points (u1 , u2 , u3 ) = (±1, 0, 0) and they will be studied in the next section. The condition b c < 0 implies that κ ∈ (2π/3, 4π/3) ∪ (−π/3, π/3), which means that e1 is the smallest or the greatest of the three inertia parameters. Dynamics 2023, 3 612 The differential Equation (32) is separated and can be integrated explicitly. Once we have r (ψ), we can substitute it into Equation (33) and solve it explicitly, or we use the first integral (31) to obtain u1 (ψ). Integrating (32) we have s 2H r (ψ) = , (35) c − b + (b + c) cos(2ψ) and for u1 q u1 ( ψ ) = ± 1 − r 2 ( ψ ), (36) where 2H, the other first integral, depends on the initial conditions. Factor 2 is convenient as will be seen ahead. The solutions r (ψ) and u1 (ψ) are well defined if (34) is satisfied, and they are periodic functions with period π. The factor 2 in the argument of the cosine function is due to the invariance of the solution under ψ → ψ + π, (u1 , u2 , u3 ) → (u1 , −u2 , −u3 ). In other words, the two points (u1 , u2 , u3 ) and (u1 , −u2 , −u3 ) are on the same solution. The first integral 2H can be easily calculated, and it must be consistent with (20). From (35) and (29), we have 2H = (c − b)(u22 + u23 ) + (b + c)r2 (cos2 ψ − sin2 ψ) = 2 c u22 − 2 b u23 , and using the definitions of b and c, (27) and (28), and the first integral (31) we obtain H = e1 u21 + e2 u22 + e3 u23 − e1 (u21 + u22 + u23 ) = e1 u21 + e2 u22 + e3 u23 − e1 = const., then e1 u21 + e2 u22 + e3 u23 = e0 = const., which is the hyperboloid of energy (20). Then, H and the energy parameter e0 are related by H = e0 − e1 . (37) We want to relate the first integral H with the values of r (0) and r (π/2). If we denote r0 = r (0) and rπ/2 = r (π/2), from (35) with ψ = 0 we have H = c r02 = −b rπ/2 2 , In this way, we have that the distances to the cylindrical axis along u2 and u3 are e0 − e1 e0 − e1 r02 = , 2 rπ/2 = , (38) e2 − e1 e3 − e1 respectively, where 0 < r02 < 1 and 0 < rπ/2 2 < 1. Thus, we have the relations e0 − e1 e0 − e1 0< < 1, 0< < 1. (39) e2 − e1 e3 − e1 From (34), we consider the following two cases. (i) If b < 0 and c > 0, then κ ∈ (2π/3, 4π/3) and e1 is the smallest inertia parameter. The denominators of (39) are positive, and for e0 we obtain e1 < e0 , e0 < e2 , e0 < e3 , (40) so that e0 is less than both e2 and e3 in all the interval κ ∈ (2π/3, 4π/3), and if (a) a < 0 then κ ∈ (2π/3, π ), e1 < e0 < e3 < e2 . Dynamics 2023, 3 613 (b) a > 0 then κ ∈ (π, 4π/3), e1 < e0 < e2 < e3 . (c) a = 0 then κ = π and the body is symmetric such that 1 − 1 = e1 < e0 < e2 = e3 = . 2 (ii) If b > 0 and c < 0, then κ ∈ (−π/3, π/3) and e1 is the greatest inertia parameter. The denominators of (39) are negative, and for e0 we obtain e1 > e0 , e0 > e2 , e0 > e3 , (41) so that e0 is greater than both e2 and e3 in all the interval κ ∈ (−π/3, π/3), and if (a) a > 0 then κ ∈ (−π/3, 0), e2 < e3 < e0 < e1 . (b) a < 0 then κ ∈ (0, π/3), e3 < e2 < e0 < e1 . (c) a = 0 then κ = 0 and the body is symmetric satisfying 1 − = e3 = e2 < e0 < e1 = 1. 2 These cases and conditions on e1 , e2 , e2 , and e0 , are shown in the Figure 1 as 1ia, 1ib, 1ic and 1iia, 1iib, 1iic. As c − b = −3e1 and b + c = e2 − e3 , the solutions (35) and (36) in terms of the inertia and energy parameters are s 2 ( e0 − e1 ) r (ψ) = , −3e1 + (e2 − e3 ) cos(2ψ) (42) p u1 ( ψ ) = ± 1 − r 2 ( ψ ), where e0 satisfies (40) or (41). As a function of the asymmetry parameter κ, Equation (13) yields s 2(e0 − cos κ )) r (ψ) = √ , (43) −3 cos κ + 3 sin κ cos(2ψ) and q u1 ( ψ ) = ± 1 − r 2 ( ψ ), which only depend on the two constant parameters κ and e0 . In the special symmetric case (1i c) with κ = π ∈ Ks , the solutions are reduced to r r 2 1 r= ( 1 + e0 ) , u1 = ± (1 − 2e0 ), (κ = π ) 3 3 which are two circles perpendicular to the u1 axis with the center at (u1 , 0, 0) and radio r. In the case (1ii c) with κ = 0 ∈ Ks , the corresponding values are r r 2 1 r= ( 1 − e0 ) , u1 = ± (1 + 2e0 ), (κ = 0). 3 3 Finally, from (38) and (29) at ψ = 0 and ψ = π/2, the trajectory crosses the planes u3 = 0 and u2 = 0, at the points e2 − e0 e0 − e1  r r  ( u 1 , u 2 , u 3 ) ψ =0 = ± , ,0 ,  r e2 − e1 e2 − e1  e3 − e0 e0 − e1 r (u1 , u2 , u3 )ψ=π/2 = ± , 0, , e3 − e1 e3 − e1 Dynamics 2023, 3 614 respectively, such that either (40) or (41) is satisfied. Case 2. Cylindric coordinates (u2 , r, ψ) around the u2 axis u1 = r sin ψ, u2 = u2 , u3 = r cos ψ. (44) The Euler Equation (25) becomes u1 u˙1 + u3 u˙3 ṙ = = ( a + c) r u2 sin ψ cos ψ, r u3 u˙1 − u1 u˙3 = u2 a cos2 ψ − c sin2 ψ ,  ψ̇ = 2 r u˙2 = b r2 sin ψ cos ψ. Now r ṙ + u2 u̇2 = 0 gives the first integral u22 + r2 = const., where r2 = u21 + u23 , which is the conservation of the norm of the angular momentum (21), i.e., u22 + r2 = 1. (45) The differential system with the independent variable ψ becomes r (ψ)( a + c) sin ψ cos ψ r 0 (ψ) = , (46) a cos2 ψ − c sin2 ψ b r2 (ψ) sin ψ cos ψ u20 (ψ) = . (47) u2 (ψ) a cos2 ψ − c sin2 ψ In order that these equations be well defined, the denominators of (46) and (47) must be different from zero for all ψ. We assume that u2 ( ψ ) 6 = 0 , a c < 0, (48) and also r < 1 by (45). When r = 0, there are two critical points (u1 , u2 , u3 ) = (0, ±1, 0); they will be studied in the next section. The condition a c < 0 implies κ ∈ (π/3, π ) ∪ (4π/3, 2π ), thus e2 is the smallest or the greatest of the three inertia parameters. The differential Equation (46) can be integrated explicitly, and with the first integral (45) we obtain u2 (ψ) s 2H r (ψ) = , (49) a − c + ( a + c) cos(2ψ) q u2 ( ψ ) = ± 1 − r 2 ( ψ ), (50) the other first integral being H. The solution is invariant under ψ → ψ + π, so the points (u1 , u2 , u3 ) and (−u1 , u2 , −u3 ) belong to the same solution. Following the same steps as in Case 1, we obtain analogous results with the corre- sponding permuted indices. For H, we obtain H = e0 − e2 and denoting r0 = r (0) and rπ/2 = r (π/2), we have H = a r02 = −c rπ/2 2 , Dynamics 2023, 3 615 therefore the square distances to the cylindrical axis u2 along u1 and u3 are given by e0 − e2 e0 − e2 r02 = , 2 rπ/2 = , (51) e3 − e2 e1 − e2 respectively, where 0 < r02 < 1 and 0 < rπ/2 2 < 1. Then, e0 − e2 e0 − e2 0< < 1, 0< < 1. (52) e3 − e2 e1 − e2 We have two cases (i) If a < 0 and c > 0, then κ ∈ (π/3, π ), and e2 is the greatest inertia parameter for all this interval of κ. The denominators in (52) are negative, and for e0 we obtain e2 > e0 , e0 > e1 , e0 > e3 , (53) so that e0 is greater than both e1 and e3 in all the interval κ ∈ (π/3, π ), and if (a) b > 0 then κ ∈ (π/3, 2π/3), e3 < e1 < e0 < e2 . (b) b < 0 then κ ∈ (2π/3, π ), e1 < e3 < e0 < e2 . (c) b = 0 then κ = 2π/3 and the body is symmetric and 1 − = e3 = e1 < e0 < e2 = 1. 2 (ii) If a > 0 and c < 0, then κ ∈ (4π/3, 2π ), and e2 is the smallest inertia parameter in all this interval of κ. The denominators of (52) are positive, and for e0 we obtain e2 < e0 , e0 < e1 , e0 < e3 , (54) so that e0 is smaller than both e1 and e3 in all the interval κ ∈ (4π/3, 2π ), and if (a) b < 0 then κ ∈ (4π/3, 5π/3), e2 < e0 < e1 < e3 . (b) b > 0 then κ ∈ (5π/3, 2π ), e2 < e0 < e3 < e1 . (c) b = 0 then κ = 5π/3 and the body is symmetric and 1 − 1 = e2 < e0 < e1 = e3 = . 2 These cases and conditions on e1 , e2 , e2 , and e0 , are shown in the Figure 1 as 2ia, 2ib, 2ic and 2iia, 2iib, and 2iic. As a − c = −3e2 and a + c = e3 − e1 , the solutions (49) and (50) in terms of the inertia and energy parameters are s 2 ( e0 − e2 ) r (ψ) = , −3e + (e3 − e1 ) cos(2ψ) (55) p 2 u2 ( ψ ) = ± 1 − r 2 ( ψ ), such that e0 satisfies (53) or (54). With (13), the previous solutions can be written in terms of only two parameters κ and e0 . In the two symmetrical cases (2i c) and (2ii c), the solutions (55) and are reduced, respectively, to r r 2 1 r= ( 1 − e0 ) , u2 = ± (1 + 2e0 ) , (κ = 2π/3) 3 3 and r r 2 1 r= ( 1 + e0 ) , u2 = ± (1 − 2e0 ), (κ = 5π/3) 3 3 Dynamics 2023, 3 616 which are circles perpendicular to the u2 axis with center at (0, u2 , 0) and radio r. Finally, from (51) and (44) at ψ = 0 and ψ = π/2, the solution crosses the planes u1 = 0 and u3 = 0 at the points e3 − e0 e0 − e2  r r  ( u 1 , u 2 , u 3 ) ψ =0 = 0, ± , , e3 − e2r e3 − e2  e0 − e2 e1 − e0 r (u1 , u2 , u3 )ψ=π/2 = ,± ,0 , e1 − e2 e1 − e2 respectively, and either (53) or (54) is satisfied. Case 3. Cylindric coordinates around the u3 axis u1 = r cos ψ, u2 = r sin ψ, u3 = u3 . (56) The Euler Equation (25) becomes u1 u˙1 + u2 u˙2 ṙ = = ( a + b) r u3 sin ψ cos ψ, r u1 u˙2 − u2 u˙1 = u3 b cos2 ψ − a sin2 ψ ,  ψ̇ = 2 r u˙3 = c r2 sin ψ cos ψ. Since r ṙ + u3 u̇3 = 0, we have the integral of the norm u23 + r2 = 1, (57) and the differential system with the independent variable ψ is r (ψ)( a + b) sin ψ cos ψ r 0 (ψ) = , (58) b cos2 ψ − a sin2 ψ c r2 (ψ) sin ψ cos ψ u30 (ψ) = . u3 (ψ) b cos2 ψ − a sin2 ψ Analogously to the previous cases, we assume that u3 (ψ) 6= 0, a b < 0, (59) and r < 1 by (57). The two critical points take place when r = 0 at (u1 , u2 , u3 ) = (0, 0, ±1). They will be studied in the next section. The condition a b < 0 means κ ∈ (0, 2π/3) ∪ (π, 5π/3) , which implies that e3 is the smallest or the greatest of the three inertia parameters. The Equation (58) is separable, and we obtain s 2H r (ψ) = , (60) b − a + ( a + b) cos(2ψ) and q u3 ( ψ ) = ± 1 − r 2 ( ψ ). (61) Both are well defined if (59) is satisfied. For H in this case, we obtain H = e0 − e3 , Dynamics 2023, 3 617 If r0 = r (0) and rπ/2 = r (π/2), we have H = br02 = − arπ/2 2 , then r r H H r0 = , rπ/2 = , b −a with 0 < r02 < 1 and 0 < rπ/2 2 < 1. Therefore e0 − e3 e0 − e3 r02 = , 2 rπ/2 = , (62) e1 − e3 e2 − e3 and e0 − e3 e0 − e3 0< < 1, 0< < 1. (63) e1 − e3 e2 − e3 We have the following two cases: (i) If a < 0 and b > 0, then κ ∈ (0, 2π/3) and e3 is the smallest inertia parameter for all this interval of κ. The denominators in (63) are positive, and for e0 we obtain e3 < e0 e0 < e1 , e0 < e2 , (64) so that e0 is smaller than both e1 and e2 in all the interval κ ∈ (0, 2π/3), and if (a) c < 0 then κ ∈ (0, π/3), e3 < e0 < e2 < e1 . (b) c > 0 then κ ∈ (π/3, 2π/3), e3 < e0 < e1 < e2 . (c) c = 0 then κ = π/3. The rigid body is symmetric with 1 − 1 = e3 < e0 < e1 = e2 = . 2 (ii) If a > 0 and b < 0, then κ ∈ (π, 5π/3), and e3 is the greatest inertia parameter for all this interval of κ. The denominators in (63) are negative, and for e0 we have e0 < e3 e1 < e0 , e2 < e0 , (65) so that e0 is greater than both e1 and e2 in all the interval κ ∈ (π, 5π/3), and if (a) c > 0 then κ ∈ (π, 4π/3), e1 < e2 < e0 < e3 . (b) c < 0 then κ ∈ (4π/3, 5π/3), e2 < e1 < e0 < e3 . (c) c = 0 then κ = 4π/3 and the rigid body is symmetric with 1 − = e1 = e2 < e0 < e3 = 1. 2 These cases and conditions on e1 , e2 , e2 , and e0 , are shown in the Figure 1 as 3ia, 3ib, 3ic and 3iia, 3iib, 3iic. As b − a = −3e3 and a + b = e1 − e2 , the solutions (60) and (61) become s 2 ( e0 − e3 ) r (ψ) = , −3e + (e1 − e2 ) cos(2ψ) (66) p 3 u3 ( ψ ) = ± 1 − r 2 ( ψ ), where e0 satisfies (64) or (65). With (13), the previous equations can be written in terms of κ, e0 and of course ψ. In the two special symmetrical cases (3i c) and (3ii c), the solutions (66) are reduced to r r 2 1 r= ( 1 + e0 ) , u3 = ± (1 − 2e0 ), (κ = π/3) 3 3 Dynamics 2023, 3 618 and r r 2 1 r= ( 1 − e0 ) , u3 = ± (1 + 2e0 ), (κ = 4π/3) 3 3 respectively. They are circles perpendicular to the u3 axis with the center at (0, 0, u3 ) and radio r. Finally, from (62) and (56) at ψ = 0 and ψ = π/2, the trajectory crosses the planes u2 = 0 and u1 = 0 at the points e0 − e3 e1 − e0 r r  (u10 , u20 , u30 ) = , 0, ± , e1 − e3 e1 − e3 e0 − e3 e2 − e0  r r  (u1π/2 , u2π/2 , u3π/2 ) = 0, ,± , e2 − e3 e2 − e3 respectively, and either (64) or (65) is satisfied. We have seen that each case has two open disjointed intervals of κ (subcases (i) and (ii)), and that the parameter space is divided in six regions where the periodic solutions were found, depending on the value of e0 . In Figure 1, we summarize the studied cases. A file Free Rigid Body can be downloaded where the solutions at the intersection of the sphere and hyperboloid can be seen as a function of κ and e0 ; it requires free CDF Wolfram Player. The boundary of these regions are the special cases where e0 = ei for some i = 1, 2, 3. We study these special cases in the next section. It must be clear from the previous discussion that the axis of the hyperboloid E , which changes at the values κ ∈ Ka (19), may not be the same as the axis of the cylindrical coordinates used to write the solution, summarized in Figure 1. For instance, for negative e0 = −0.5 and κ ∈ (2π/3, 4π/3) (cases (1ia–c)), when κ ∈ (2π/3, 5π/6), E is a one sheet hyperboloid with axis u2 , but the solution has the cylindrical axis u1 (case (1ia)). For κ ∈ (5π/6, 7π/6), E is a two sheet hyperboloid with axis u1 , and the solution has the same cylindrical axis u1 (cases (1ia–c)). For κ ∈ (7π/6, 4π/3), E is a one sheet hyperboloid with axis u3 , but the solution has the cylindrical axis u1 (case (1ib)). Finally, with the division into six disjoint regions based on the symmetry axis of the cylindrical coordinates to obtain the solutions, we can see that the values of e0 and κ are related to the solutions found in each case. The cylindrical symmetry axes where the solutions were calculated coincide with the axes of the larger or smaller inertia parameters, never with the intermediate one. For each value of κ, there is always an intermediate inertia parameter, which we will refer to as eint , and two extremes, denoted as emin < 0 < emax . If the energy parameter e0 satisfies emin < e0 < eint < emax , then the solution is expressed in cylindrical coordinates with the cylindrical axis in the direction of emin . Conversely, if the energy parameter e0 satisfies emin < eint < e0 < emax , then the solution is expressed in cylindrical coordinates with the cylindrical axis in the direction of emax . In the symmetric cases, there are two equal inertia parameters, let us say ei = e j . In this scenario, the different inertia parameter ek satisfies either ek < e0 < ei = e j or ei = e j < e0 < ek . These symmetric cases allow the exchange of the cylindrical axis in the i and j directions for the calculated solution: a parameter ei that was intermediate for κ slightly lower than the value of the symmetric case (ei < e j ), becomes extreme for κ slightly greater than in the symmetric case (e j < ei ), and the other one that was extreme becomes intermediate. Of course, the third parameter ek remains extreme. However, the axis of the cylindrical coordinates for the solution may not necessarily correspond with the axis of the hyperboloid, as we mentioned earlier. 4. The Separatrices Now we consider the cases in which e0 takes the value of an inertia parameter ei at the boundary of the regions studied in the previous section. The most simple cases are those in which r = 0 because the solution reduces to two points. This is confirmed by calculating the limit of the solutions previously found, (42), (55) and (66), when the numerator goes to Dynamics 2023, 3 619 zero. For instance, in Case 1, we saw that e1 has the smallest or the greatest value of the three inertia parameters, and from (42), if e0 = e1 , then r = 0 and (u1 , u2 , u3 ) = (±1, 0, 0). These are two points on the u1 axis, where the two surfaces (20) and (21) are tangent, as it is expected. However, in the limit when e0 tends to the inertia parameter with the intermediate value (for a given κ), the solutions are the separatrix sets given by four semicircles (they are open sets). To see this we substitute the value of e0 in (42), (55) and (66) (Cases 1, 2 and 3, respectively), to the intermediate value of the inertia parameter, as it is indicated in the second column of the next tables. √ Case 1 u1 = ± 1 − r 2 (i a) κ ∈ (2π/3, π ) r √   −3 cos √ κ − 3 sin κ  π 3π e0 = e3 , r = , ψ ∈ − π2 , π 2 ∪ 2, 2 (ii a) κ ∈ (5π/3, 2π ) −3 cos κ + 3 sin κ cos 2ψ (i b) κ ∈ (π, 4π/3) r √ −3 cos √ κ + 3 sin κ e0 = e2 , r = , ψ ∈ (0, π ) ∪ (π, 2π ) (ii b) κ ∈ (0, π/3) −3 cos κ + 3 sin κ cos 2ψ (i c) κ = π e0 = e2 = e3 , r=1 (ii c) κ = 0 √ Case 2 u2 = ± 1 − r 2 r (i a) κ ∈ (π/3, 2π/3) 2(cos(κ )−cos(κ −2π/3)) e0 = e1 , r = −3 cos(κ −2π/3)−(cos κ −cos(κ +2π/3)) cos(2ψ) ,    (ii a) κ ∈ (4π/3, 5π/3) ψ ∈ − π2 , π2 ∪ π2 , 3π 2 r (i b) κ ∈ (2π/3, π ) 2(cos(κ +2π/3)−cos(κ −2π/3)) e0 = e3 , r = −3 cos(κ −2π/3)−(cos κ −cos(κ +2π/3)) cos(2ψ) , (ii b) κ ∈ (5π/3, 2π ) ψ ∈ (0, π ) ∪ (π, 2π ) (i c) κ = 2π/3 e0 = e1 = e3 , r=1 (ii c) κ = 5π/3 √ Case 3 u3 = ± 1 − r 2 r (i a) κ ∈ (0, π/3) 2(cos(κ −2π/3)−cos(κ +2π/3)) e0 = e2 , r = −3 cos(κ +2π/3)+(cos(κ )−cos(κ −2π/3)) cos(2ψ) ,   π 3π  (ii a) κ ∈ (π, 4π/3) ψ ∈ − π2 , π 2 ∪ 2, 2 r (i b) κ ∈ (π/3, 2π/3) 2(cos(κ )−cos(κ +2π/3)) e0 = e1 , r = −3 cos(κ +2π/3)+(cos(κ )−cos(κ −2π/3)) cos(2ψ) , (ii b) κ ∈ (4π/3, 5π/3) ψ ∈ (0, π ) ∪ (π, 2π ) (i c) κ = π/3 e0 = e1 = e2 , r=1 (ii c) κ = 4π/3 In general, assuming that the intermediate value of the inertia parameters is e j , and their ordering is emin = ei < e j < ek = emax , i 6= j 6= k in the set {1, 2, 3}, if e0 = e j , the solutions are four open semicircles whose boundary points are on opposite sides of the u j axis. The intersection of the hyperboloid (20) and the unit sphere (21) are two big circles, where the four semicircles lie. Since the separatrix solution is at the boundary of the two different cases, it can be obtained indistinctly from either of them, as the locus of both solutions coincides, although the domain of ψ is different. For instance, solutions (1ia) are the same locus than (2ib), and (1iia) are the same locus than (2iib), although they come from different limits and symmetry axes. The domain of (1ia) and (1iia) is ψ ∈ − π2 , π2 ∪ π2 , 3π  2 and of the second one is ψ ∈ (0, π ) ∪ (π, 2π ). 5. Conclusions This work was intended to provide an entirely different approach to the torque free Euler equations by mean of the reduction to a cylindrical variable ψ as the new independent Dynamics 2023, 3 620 variable of the differential system. We used the dimensionless variables and parameters to characterize the rigid body with only two parameters instead of five (three moments of inertia, E and L); they are the asymmetry parameter κ and the energy-momentum parameter e0 . We made the reduction of the differential system eliminating the time as the independent variable in the Euler equations in dimensionless variables, and we studied the solution and its geometrical aspects by means of continuous deformations dependent on the two independent parameters κ and e0 . The parameter space was divided into six disjoint regions, whose boundaries are the limit cases where the symmetry axis changes. The solutions depend on where in the six regions—or their boundaries—the values of e0 and κ lie in the parameter space, as in each of them an specific symmetry axis is convenient to write the solution. The cases named 1 have the cylindrical axis u1 , the cases 2 have the cylindrical axis u2 , and the cases 3 have the cylindrical axis u3 . Continuously changing the two independent parameters, we obtained the solution r (ψ) and ui (ψ) in cylindrical coordinates, with ui the cylindrical axis appropriate in each case. The solutions at the boundaries of the regions are the separatrices where there are two symmetry axis, which allows the transition from a symmetry axis to another. Four open circles were found explicitly as a function of κ, such that e0 takes the value of the intermediate inertia parameter. Author Contributions: Both authors contributed equally to this paper. All authors have read and agreed to the published version of the manuscript. Funding: The second author is partially supported by the Ministerio de Economía, Industria y Competitividad, Agencia Estatal de Investigación grant MTM 2016-77278-P (FEDER) and PID2019- 104658GB-I00 (FEDER), the Agència de Gestió d’Ajuts Universitaris i de Recerca grant 2017 SGR 1617, and the European project Dynamics-H2020-MSCA-RISE-2017-777911. Data Availability Statement: Not applicable. Conflicts of Interest: The authors declare no conflict of interest. References 1. Marquina, J.; Marquina, M.; Marquina, V.; Hernández-Gómez, J. Leonhard Euler and the mechanics of rigid bodies. Eur. J. Phys. 2017, 38, 015001. [CrossRef] 2. Poinsot, L. Théorie Nouvelle de la Rotation des Corps; Libraire de Bachelier: Paris, France, 1851. 3. Rueb, A.S. Specimen Inaugurale de Motu Gyratorio Corporis Rigidi, Nulla vi Acceleratrice Sollicitati. Ph.D. Thesis, Academia Rehno-Trajectina, Utrecht, The Netherlands, 1834. 4. Jacobi, C.G. Sur la rotation d’un corps. Extrait d’une lettre adressée à l’académie des sciences de Paris. J. Die Reine Angew. Math. 1850, 39, 293–350. 5. Whittaker, E.; Watson, G. A Course of Modern Analysis: An Introduction to the General Theory of Infinite Processes and of Analytic Functions; with an Acount of the Principal Transcendental Functions, 4th ed.; Cambridge University Press: Cambridge, UK, 1996. 6. Whittaker, E. A Treatise on the Analytical Dynamics of Particles and Rigid Bodies, 4th ed.; Cambridge University Press: Cambridge, UK, 1965. 7. Lagrange, J.L. Mechanique Analitique; La Veuve Desaint: Paris, France, 1788. 8. Kovalevskaya, S. Sur le problème de la rotation d’un corps solide autour d’un point fixe. Acta Math. 1889, 12, 77–232. 9. Goryachev, D. On a motion of a heavy rigid body about a fixed point in the case of A = B = 4C. Mat. Sbonik Kruzhka Lyub. Mat. Nauk 1900, 21, 431–438. 10. Leimanis, E. The General Problem of the Motion of Coupled Rigid Bodies about a Fixed Point; Springer Science & Business Media: Berlin/Heidelberg, Germany, 2013; Volume 7. 11. Bechlivanidis, C.; Van Moerbeke, P. The Goryachev-Chaplygin top and the Toda lattice. Commun. Math. Phys. 1987, 110, 317–324. [CrossRef] 12. Arnold, V.I. A theorem of Liouville concerning integrable problems of dynamics. Sib. Mat. Zhurnal 1963, 4, 471–474. 13. Arnold, V.I. Mathematical Methods of Classical Mechanics; Springer Science & Business Media: Berlin/Heidelberg, Germany, 2013; Volume 60. 14. Andoyer, H. Mécanique Céleste; Gauthier-Villars: Paris, France, 1926. 15. Deprit, A. Free rotation of a rigid body studied in the phase plane. Am. J. Phys. 1967, 35, 424–428. [CrossRef] 16. Marsden, J.; Ratiu, T. Introduction to Mechanics and Symmetry, 2nd ed.; Springer: Berlin/Heidelberg, Germany, 2003. 17. Holm, D.D. Geometric Mechanics-Part I: Dynamics and Symmetry; Imperial College Press: London, UK, 2011. 18. Holm, D.D.; Marsden, J.E. The rotor and the pendulum. In Symplectic Geometry and Mathematical Physics; Birkhäuser: Boston, MA, USA, 1991; pp.189–203. Dynamics 2023, 3 621 19. de la Cruz, M.; Gaspar, N.; Jiménez-Lara, L.; Linares, R. Classification of the classical SL (2, R) gauge transformations in the rigid body. Ann. Phys. 2017, 379, 112–130. [CrossRef] 20. Reiche, F. Zur Quantelung des asymmetrichen Kreisels. Phys. Z. 1918, 19, 394–399. 21. King, G.W. The asymmetric rotor. VI. Calculation of higher energy levels by means of the correspondence principle. J. Chem. Phys. 1947, 15, 820–830. [CrossRef] 22. Kramers, H.; Ittmann, G. Zur quantelung des asymmetrischen kreisels. Z. Phys. 1929, 53, 553–565. [CrossRef] 23. Lukac, I.; Smorodinskiı̌, Y.A. The wave functions of an asymmetric top. Sov. J. Exp. Theor. Phys. 1969, 30, 728. 24. Spence, R. Angular momentum in sphero-conal coordinates. Am. J. Phys. 1959, 27, 329–335. [CrossRef] 25. Piña, E. Solution by the separation method of the motion of a rigid body with no forces. Rev. Mex. Física 1997, 43, 205–218. 26. Ley-Koo, E. Rotations of the most asymmetric molecules via 4-step and 1-step ladder operators. J. Phys. Conf. Ser. 2015, 597, 012055. [CrossRef] 27. Borisov, A.; Mamaev, I. Rigid Body Dynamics; Higher Education & Press Walter de Gruyter GmBH: Berlin, Germany; Boston, MA, USA, 2019. 28. Chernousko, F.L.; Akulenko, L.D.; Leshchenko, D.D. Evolution of Motions of a Rigid Body about Its Center of Mass; Springer: Berlin/Heidelberg, Germany, 2017. 29. Murakami, C. Analytical Solution of the Euler-Poinsot Problem. J. Geom. Symmetry Phys. 2021, 60, 25–46. [CrossRef] 30. Pezent, J.B.; Sood, R. New Solution to Euler’s Equations for Asymmetric Bodies: Applications to Spacecraft Reorientation. J. Guid. Control. Dyn. 2022, 45, 2289–2303. [CrossRef] 31. Yehia, H.M. Rigid Body Dynamics: A Lagrangian Approach; Springer Nature: Cham, Switzerland, 2022. 32. Zhuravlev, V.F.; Rozenblat, G.M. Estimates of Solutions During Motion of the Euler-Poinsot Top and Explanation of the Experiment with Dzhanibekov’s Nut. Russ. J. Nonlinear Dyn. 2020, 16, 517–525. [CrossRef] 33. Piña, E. Dinámica de Rotaciones; Colección CBI, Universidad Autónoma Metropolitana: Mexico City, Mexico, 1996. 34. Landau, L.; Lifshitz, E. Mechanics, 3rd ed.; Pergamon Press: Oxford, UK, 1989; Volume 1. Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.