Study of Stability of Rotational Motion of Spacecraft with Canonical Variables

1 Group of Orbital Dynamics and Planetology, São Paulo State University (UNESP), Guaratinguetá 12516-410, SP, Brazil 2 Department of Mathematics, São Paulo Salesian University (UNISAL), Lorena 12600-100, SP, Brazil 3 Space Mechanic and Control Division, National Institute for Space Research (INPE), São José dos Campos 12227-010, SP, Brazil 4 Department of Aircraft Maintenance and Aeronautical Manufacturing, Faculty of Technology (FATEC), São José dos Campos 12247-004, SP, Brazil


Introduction
Stability analysis of the rotational motion of a satellite taking into account the influence of external torques is very important in maintaining the attitude to ensure the success of a space mission.
Recently, some studies on the subject have been developed, and they motivated the development of this paper.In 1, 2 is presented a study on the stability of the rotational motion of artificial satellites in an elliptic orbit, considering disturbance due to gravity gradient torque, using canonical formulations and the Andoyer canonical variables.These studies use the procedure presented in 3 to determine a normal form of the Hamiltonian up to 4th order.
In 4 is developed a numerical-analytical method for normalization of Hamiltonian systems with 2 and 3 degrees of freedom.The normal form is obtained using the Lie-Hori method 5 .The stability analysis of the system is made by the Kovalev-Savchenko theorem 6 .The most important role of this work is the results obtained analytically for the generating function of 3rd order, necessary for determining the coefficients of the normal Hamiltonian of 4th order.
Thus, the objective of this work is to optimize the stability analysis developed in 2 to determine the equilibrium points, the normal form for dynamical systems with two degrees of freedom, and the applications of the Kovalev-Savchenko theorem 6 .It will be done by applying the expressions obtained in 4 for the coefficients of the normal 4th-order Hamiltonian.
In this paper equilibrium points and/or regions of stability are established when parcels associated with gravity gradient torque acting on the satellite are included in the equations of rotational motion.
The Andoyer variables are used to describe the rotational motion of the satellite in order to facilitate the application of methods of stability of Hamiltonian systems.The Andoyer canonical variables 7 are represented by generalized moments L 1 , L 2 , L 3 and by generalized coordinates 1 , 2 , 3 that are outlined in Figure 1.The angular variables 1 , 2 , 3 are angles related to the satellite system Oxyz with axes parallel to the spacecraft's principal axes of inertia and equatorial system OXYZ with axes parallel to the axis of the Earth's equatorial system .Variables metrics L 1 , L 2 , L 3 are defined as follows: L 2 is the magnitude of the angular momentum of rotation L 2 , L 1 is the projection of L 2 on the z-axis of principal axis system of inertia L 1 L 2 cos J 2 , where J 2 is the angle between the z-satellite axis and L 2 , and L 3 is the projection of L 2 on the Z-equatorial axis L 3 L 2 cos I 2 , where I 2 is the angle between Z-equatorial axis and L 2 .
The nonlinear stability of equilibrium points of the rotational motion is analyzed here by the Kovalev-Savchenko theorem 6 , which requires the normalized Hamiltonian up to terms of fourth order around the equilibrium points.
The equilibrium points are found from the equations of motion.With the application of the Kovalev-Savchenko theorem, it is possible to verify if they remain stable under the influence of terms of higher order of the normal Hamiltonian.
In this paper numerical simulations were made for two hypothetical groups of artificial satellites, considering them in a circular orbit and with symmetric shape in relation to their physical and geometric characteristics.The satellites are classified as medium and small sized; they have orbital data and physical characteristics similar to real satellites.

Equations of Motion
The Andoyer variables, defined above, are used to characterize the rotational motion of a satellite around its center of mass 7 , and the Delaunay variables describe the translational motion of the center of mass of the satellite around the Earth 8 .
The Delaunay variables L, G, H, l, g, h are defined as 8 of the ascending node, M is the mass of the satellite, μ is the Earth gravitational parameter, I is the inclination of the orbit, a is the semimajor axis, and e is the orbital eccentricity.
Here it is assumed that the satellites are in a circular orbit, which differs from the study presented in 2 which considers satellites in eccentric orbits.This consideration was adopted to simplify the Hamiltonian of the problem, which is extensive 1 , and to facilitate the stability analysis of the equilibrium points.
Thus, assuming that the satellites have well-defined circular orbit, the main focus of the work is only aimed at stability analysis of the rotational movement of the satellites in study.
In this case the Hamiltonian of the problem is expressed in terms of the Andoyer and Delaunay variables 9 as follows: where F o is the unperturbed Hamiltonian and F 1 is the term of the Hamiltonian associated with the disturbance due to the gravity gradient torque, both are described respectively by 10 where m 2, 3 and n 1, 2, 3; A, B, and C are the principal moments of inertia of the satellite on x-axis, y-axis, and z-axis, respectively; H 1 and H 2 are functions of the variables m , L n , where 2 and 3 appear in the arguments of cosines.The complete analytical expression is presented in 1 for eccentric orbits.In this paper F 1 will be simplified considering the satellite in a circular orbit; its complete analytical expression is given in the Appendix.Thus, the equations of motion associated with the Hamiltonian 2.1 are given by These equations are used to find the possible equilibrium points of the rotational movement.
In this study, it is also considered that the satellite has two of its principal moments of inertia equal, B A. With this relationship, the variable 1 will not be present in the Hamiltonian, reducing the dynamic system to two degrees of freedom, a necessary condition for applying the stability theorem chosen for analysis of equilibrium points.

Normal Form of the Hamiltonian of 2nd Order and Linear Stability Analysis
Because the variable 1 is a cyclic coordinate, the Hamiltonian 2.1 is written as F L 2 , L 3 , 2 , 3 and the equations of motion 2.3 can be written in vector form as 1 ẇ JF w , 3.1 where w is the state vector and F w is the matrix of partial derivatives with respect to their respective variables given by and J a symplectic matrix given by Now, it is necessary to linearize the system around the equilibrium point.However, it is convenient to make a translation of the coordinates in the Hamiltonian 2.1 so that the origin coincides with the equilibrium point under study; it means that where L 1e , 2e , L 2e , 3e , and L 3e are the coordinates of equilibrium points.
With this translation we can expand the Hamiltonian in Taylor series around the new origin that is nothing more than expand the Hamiltonian in the neighborhoods of the equilibrium point when q 1 p 1 q 2 p 2 0.Then, where m 2, 3, n 1, 2, 3, and F k is the Hamiltonian expanded up to terms of order k, with k 2, 3, 4, . . . .The Hessian P of the problem is calculated from the Hamiltonian expanded up to 2nd order around the equilibrium point: Thus, the system of linearized equations can be written as where W is the state vector in the new variables So, to find the eigenvalues of the reduced linear system, we should diagonalize the matrix JP present in 3.7 by means of the determinant det λI −JP 0, where I is an identity matrix and λ are the eigenvalues of the matrix JP .If these eigenvalues are pure imaginary λ j ±iω j , the normal form of 2nd-order Hamiltonian can be expressed as 2 3.9

Extension of the Normal Form of the Hamiltonian to Higher Orders
Extending the process of normalization of the Hamiltonian to higher orders H 3 , H 4 , . . ., the Hamiltonian normal in terms of variables q i and p i can be expressed by To obtain the Hamiltonian normal H * 3 , H * 4 , . .., it is necessary to use the method of Lie-Hori.According to 11 , it is desirable that the normal Hamiltonian is written in complex variables before using this method.
This transformation of variables is given by with inverse given by

3.12
Thus, the normal Hamiltonian of 2nd order is given by iω j x j y j 3.13 and normal Hamiltonian until the 4th order in complex variables is represented by

Series Method of Lie-Hori
The normal form is obtained from the expansion of the Hamiltonian in terms of the Lie series given by 12, 13

3.22
As H new 2 is already given in its normal form, the next step is to calculate the generating function which leads to the normal form up to 4th-order terms.This procedure allows to obtain the minimum number of nonresonant monomials defined in 12 ; thus, the generating function G n is chosen to eliminate the complex variables in the total Hamiltonian 3.5 in which the monomials x j and y j have different exponents of the desired degree, leaving only the monomials that carry on the resonance of the intrinsic Hamiltonian systems.The following is a procedure to obtain the normal Hamiltonian to 4th order 1, 4 .
I It is known that Hamiltonian systems are naturally resonant in even orders, 2nd, 4th, and so forth 10, 12 , so to find the normal Hamiltonian until the 4th order, the 3rd disturbance H new 3 must first be removed from H new , it means H new 3 0.

3.23
Assuming H * 3 is briefly defined in 3.15 to eliminate terms of 3rd order, we consider the generating function G 3 x, y expressed in abbreviated form as G 3 x, y |α| |β| 3 g 3,α,β x α y β .

3.24
Then it is necessary to find G 3 x, y of the third-order homological equation 3.21 , where

3.25
After some algebraic manipulations, we find G 3 x, y given by 4
For the resonance condition it is not accepted that β − α, W / 0, which is always true for monomials of third-order Hamiltonian when considering only the natural resonance of the Hamiltonian.II This process is repeated for the polynomial of degree 4 finding G 4 x, y that will eliminate some of monomials H * 4 in 3.16 .The remaining monomials in H new 4 are called resonant monomials 3 .Then 3.22 can be reduced as 1

3.31
After some algebraic manipulations, we find G 4 x, y given by 4 where α − β, W / 0 represents an inner product as mentioned above.Thus, using 3.34 in 3.29 it is possible to find H new 4 after some algebraic manipulations.
This procedure can be repeated for any order that is required in order to find a normal form of the new Hamiltonian H new .In this paper it is considered the 4th order of the H new .

Normal Form of the Hamiltonian of 4th Order and Nonlinear Stability Analysis
As we get the normal form of 2nd-order Hamiltonian 3.9 , we can apply the method described previously to find the Hamiltonian normal until the 4th order using the Hamiltonian in complex variables 3.14 represented by 1, 4 where H * 3 x 1 , y 1 , x 2 , y 2 and H * 4 x 1 , y 1 , x 2 , y 2 are obtained from expressions

3.34
Using the method of Lie-Hori presented before, it is possible to determine the generating function G 3 x 1 , y 1 , x 2 , y 2 in order to satisfy the homological equation 3.25 , where G 3 x 1 , y 1 , x 2 , y 2 was determined analytically in 4 .This generating function is responsible for the elimination of the 3rd-order terms of the Hamiltonian 3.33 .
The generating function G 3 x 1 , y 1 , x 2 , y 2 is used to find the terms of 4th order, F 4 , and perform the separation described in 3.28 , obtaining H 4 and G 4 x 1 , y 1 , x 2 , y 2 .The coefficients of the normal form H new 4 can be calculated using the terms obtained in H 4 and the eigenvalues of H 2 .
Thus, the normal Hamiltonian in complex variables separated by degree can be expressed as follows 4 : where ω j j 1, 2 is the imaginary part of eigenvalues associated with the matrix defined by the product of a matrix symplectic of order 4 with the Hessian of the Hamiltonian expanded in Taylor series up to 2nd order around the equilibrium point; δ ij are real coefficients obtained from the combination of h k,α,β x α y β k 3, 4 .However, the Kovalev-Savchenko theorem requires that the normal Hamiltonian is in real variables, and with the applications of the transformation of variables 3.12 the normal form Hamiltonian can been expressed as follows:

Kovalev-Savchenko Theorem
To use the stability theorem, the normal Hamiltonian of the problem is necessary as it was discussed above.
Considering the normal Hamiltonian H o , an analytics function of coordinates q ν and generalized moments p ν to a fixed point P is expressed by where O 5 represents higher-order terms; ω o ν is the imaginary part of eigenvalues associated with the matrix defined by the product of a 4th-order matrix symplectic with the Hessian of the Hamiltonian expanded in Taylor series up to 2nd order around the equilibrium point; δ o ν,υ depend on the eigenvalues ω o ν and the coefficients of the Hamiltonian expanded in Taylor series of 3rd and 4th order around the equilibrium point and they are presented analytically in Formiga 4 , and The stability analysis is performed here by the Kovalev-Savchenko theorem 6 which ensures that the motion is Lyapunov stable if the following conditions are satisfied.
i The eigenvalues of the reduced linear system are pure imaginary ±iω o 1 and ±iω o 2 .ii The condition is valid for all k 1 and k 2 integers satisfying the inequality iii The determinant D o must satisfy the inequality where δ o ν,υ are the coefficients of the normal 4th-order Hamiltonian.
This theorem says that a Hamiltonian reduced in its normal form up to 4th order, in the absence of the resonance condition of the eigenvalues associated and if the condition 4.5 is satisfied, it is guaranteed the existence of tori invariant in a neighborhood small enough of equilibrium position 14-16 .

Computational Algorithm for the Normal Form of the Hamiltonian of Rotational Motion and Analysis of Stability
In order to synthetize the process to analyze the stability of equilibrium points, the logical sequence of the algorithm is now presented.

Stability Analysis of Equilibrium Points
See Figure 2.

Numerical Simulations
As already mentioned, we consider two types of satellites: medium sized MP , which has similar orbital characteristics with the American satellite PEGASUS 17 , and small sized PP , which has similar orbital characteristics of the Brazilian data collection satellite SCD-1 and SCD-2 18 .All the numerical simulations were developed using the software MATHEMATICA.
Table 1 shows a quantitative summary of the equilibrium points found and their stability, according to the criteria of Kovalev-Savchenko theorem 6 .It can be observed that, when the first condition is not satisfied, the eigenvalues associated with the matrix JP are real or are not pure imaginary.This equilibrium point is not linearly stable.
There were found totals of 60 equilibrium points for the MP satellite, 50 equilibrium points to the PP-1 satellite and 50 equilibrium points for the PP-2 satellite.Tables 2, 3 and 4, show two equilibrium points found in the simulations, one being Lyapunov stable and the other unstable, to satellites MP, PP1 and PP2, respectively.
For the Lyapunov stable equilibrium points of Tables 2, 3, and 4 respectively, Table 5 shows the values of the angles I 2 , J 2 , the rotation speed ω, and the rotation period T .These values characterize the nonexistence of singularities in these points; it means that the angles I 2 and J 2 are not null or close to zero.By Figure 1, when these angles I 2 and J 2 are null or close to zero, the Andoyer variables l 1 , l 2 , and l 3 are indeterminate, because it is difficult to determine the intersection between the involved planes in the definitions of these variables.This analysis was performed for all equilibrium points found in the simulations.Table 6 shows the same conditions for the unstable points of Tables 2, 3, and 4.

Conclusion
In this paper we presented a semianalytical stability of the rotational motion of artificial satellites, considering the influence of gravity gradient torque for symmetric satellite in a circular orbit.Applications in numerical simulations, performed with the software MATHEMATICA, were made to two types of satellites: medium MP and small PP .
Initially the points of equilibrium were determined using the physical, orbital, and attitude characteristics of each satellite.Then the algorithm for stability analysis was applied and it was obtained 2 stable equilibrium points for the satellite MP, 7 stable points for the satellite PP-1, and 10 stable points for the satellite PP-2.For the satellite MP it was gotten only two equilibrium points because this satellite has similar characteristics to the satellite PEGASUS, which is tumbling 17 .For satellites PP-1 and PP-2 were obtained many other equilibrium points, but most of them were discarded, because they lead to the Andoyer variables, a condition of uniqueness it means that the angles I 2 and J 2 are null or close to zero, and the Andoyer variables l 1 , l 2 , and l 3 are indeterminate .
It can be observed that a larger number of stable equilibrium points were determined in comparison with the results of 1, 2 , which show the stability analysis of the rotational motion with the gravity gradient torque, but with the satellite in an elliptic orbit.
An optimization was done in the algorithm of determining the normal form of the Hamiltonian and the stability analysis, using expressions obtained by 4 for the coefficients of the normal 4th-order Hamiltonian.The introduction of these expressions enabled a more effective stability analysis for equilibrium points in comparison with the results of 1, 2 .It is possible to say that the numerical simulations have become less laborious allowing analysis of data in more numbers .
This paper presents results that can directly contribute in maintaining the attitude of artificial satellites.Once the regions of stability are known for the rotational motion, a smaller number of maneuvers to maintain the desired attitude can be accomplished.In this case, a fuel economy can be generated to the satellite with propulsion system, increasing the satellite's lifetime.

Appendix
The disturbance Hamiltonian F 1 , due to the gravity gradient torque, for a satellite in a circular orbit with two of its principal moments of inertia equal, B A, is given by 19 where: A B and C are the principal moments of inertia on x-axis, y-axis and z-axis respectively, M is the mass of the satellite, μ is the gravitational constant, I is the inclination of the orbit, I 2 is the inclination of the plane of angular momentum with the plane of the equator, and J 2 is the inclination of the principal plane with the plane of the angular momentum.
In A.1 the generalized moments L 1 , L 2 , L 3 are implicit in some terms, using the definitions of generalized moments present in the introduction and trigonometric properties, they can be explicit replacing. A.2 H G cos I, l is the mean anomaly, g is the argument of perigee, h is the longitude

Figure 2 :
Figure 2: Flowchart representative of the stability analysis of equilibrium points.
16 3.15 and 3.16, is always assumed it α, β > 0. For 3.15 , |α| |β| 3 means that the sum of the exponents of the new variables must obey the order of the Hamiltonian H * 3 , and, for 3.16 , |α| |β| 4 means that the sum of the exponents of the new variables must obey the order of the Hamiltonian H * 4 .These equations are presented entirely in 1, 4 .Thus, we have the Hamiltonian written in complex variables.It is important to note that H 2 is already normalized, and now we can apply the method of Lie-Hori to find the normal form to higher orders H * 3 , H * 4 .

Table 1 :
Quantitative summary of the classification of equilibrium points.

Table 5 :
Analysis of possible singularities of the Anloyer variables for satellites MP, PP1 of and PP2.Lyapunov stable equilibrium point of Table2, 3 and 4, respectively.

Table 6 :
Analysis of possible singularities of the Anloyer variables for satellites MP, PP1 and PP2.